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ABSTRACT 


Title  of  Dissertation:  NOISY  PRECURSORS  FOR  NONLINEAR 

SYSTEM  INSTABILITY  WITH  APPLICATION 
TO  AXIAL  FLOW  COMPRESSORS 

Taihyun  Kim,  Doctor  of  Philosophy,  1997 

Dissertation  directed  by:  Professor  Eyad  H.  Abed 

Department  of  Electrical  Engineering 


This  dissertation  addresses  monitoring  of  nonlinear  systems  for  detection  and 
prediction  of  incipient  instabilities.  The  analysis  and  design  presented  here  rely 
on  the  influence  of  noise  on  system  behavior  near  the  onset  of  instability.  The 
work  is  of  relevance  to  high  performance  engineering  systems,  which  are  often 
operated  with  a  low  stability  margin  in  order  to  maximize  performance.  In  such 
a  stressed  operating  mode,  a  small  or  moderate  disturbance  can  result  in  loss  of 
stability  of  the  nominal  operating  condition.  This  can  be  followed  by  operation  in 
a  new  lower  performance  mode,  oscillatory  behavior,  or  even  system  collapse.  All 
of  these  conditions  can  be  viewed  as  bifurcations  in  the  underlying  dynamical 
models.  Prediction  of  the  precise  onset  points  of  these  instabilities  is  made 
difficult  by  the  lack  of  accurate  models  for  complex  engineering  systems.  Thus, 


in  this  thesis  monitoring  systems  are  proposed  that  can  signal  an  approaching 
instability  before  it  occurs,  without  requiring  a  precise  system  model. 

The  approach  taken  in  this  work  is  based  on  precursors  to  instability  that 
are  features  of  the  power  spectral  density  of  a  measured  output  signal.  The 
noise  in  the  system  can  be  naturally  occurring  noise  or  can  be  intentionally 
injected  noise.  The  output  signal  can  be  measured  directly  from  the  physical 
system  or  from  the  system  with  an  augmented  monitoring  system.  Design  of 
appropriate  augmented  monitoring  systems  is  a  major  topic  of  this  work.  These 
monitoring  systems  result  in  enhancing  precursor  signals  and  also  allow  control 
of  the  precursor  by  tuning  external  parameters.  This  tuning  is  important  in  that 
it  adds  confidence  to  the  detection  of  an  impending  instability. 

The  methods  developed  on  precursors  for  instablilitv  are  applied  to  models  of 
axial  flow  compression  systems.  Existing  results  on  bifurcations  for  such  models 
and  their  relation  to  compressor  stall  provide  a  starting  point  for  the  analysis. 
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Chapter  1 


Introduction 


In  recent  years,  bifurcation  and  chaos  have  been  a  very  active  research  area 
in  which  a  good  deal  of  theoretical  and  applied  work  has  been  done.  Bifur¬ 
cation  theory  has  an  important  role  in  science,  engineering  and  social  science. 
The  appearance  of  topologically  nonequivalent  phase  portraits  under  variation 
of  parameters  is  called  a  bifurcation  [38] .  Systems  that  upon  analysis  are  found 
to  be  nonlinear,  nonequilibrium,  deterministic,  dynamic,  and  that  incorporate 
randomness  so  that  they  are  sensitive  to  initial  conditions,  and  have  strange 
attractors  are  said  to  be  chaotic  [16]. 

Due  to  today’s  ever-increasing  demands  for  performance,  engineering  sys¬ 
tems  are  often  required  to  operate  very  close  to  the  limits  of  their  operating- 
envelopes.  Exceeding  these  limits  can  result  in  degradation  of  performance,  os¬ 
cillatory  behavior,  a  jump  to  an  unacceptable  new  operating  condition,  or  even 
system  collapse. 

Operation  near  the  limits  of  a  system  envelope  ( “ stressed  operation1 ' )  can 
still  lead  to  these  undesirable  effects.  This  is  because  even  a  small  or  moderate 
external  disturbance  can  push  the  system  outside  its  operating  envelope.  Air¬ 
craft  stall  in  supermaneuvers,  blackout  of  a  heavily  loaded  power  system,  and 
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compressor  stall  are  just  a  few  examples  of  system  failures  that  can  occur  in 
stressed  operation. 

Conventional  linear  models  are  incapable  of  explaining  such  catastrophic  dy¬ 
namical  behavior.  Nonlinear  dynamic  models  can  often  be  employed  to  explain 
transitions  in  system  behavior.  In  particular,  bifurcation  theory  has  proved  to 
be  quite  a  useful  tool  for  understanding  qualitative  changes  in  the  behavior  of 
nonlinear  systems.  This  theory  has  become  an  indispensable  tool  in  the  analysis 
of  many  problems  in  many  technological  and  scientific  problems.  In  [7],  [10], 
[41],  [39],  and  [43],  bifurcation  theory  has  been  successfully  applied  to  axial  flow 
compressor  system  models.  It  has  also  been  employed  to  explain  the  dynamics 
of  voltage  collapse  in  electric  power  systems  in  [23],  [4],  [12],  [18],  [17],  and  [55]. 

The  fact  that  bifurcations  are  usually  considered  as  undesirable  events  in 
applications  has  motivated  research  on  the  control  of  bifurcations.  There  are 
three  basic  categories  of  control  of  bifurcations. 

One  is  the  use  of  a  state  feedback  to  relocate  (delay)  a  bifurcation  in  pa¬ 
rameter  space.  Linear  state  feedback  control  is  often  enough  to  achieve  this. 
Application  of  linear  state  feedback  to  delay  bifurcation  in  axial  flow  compres¬ 
sors  is  given  in  [50],  [51]. 

Another  category  of  bifurcation  control  involves  modifying  the  stability  of 
bifurcated  solutions.  In  [5],  Abed  and  Fu  devised  local  feedback  laws  that  ensure 
that  the  periodic  orbit  emerging  at  a  Hopf  bifurcation  is  stable.  In  [6],  they 
extended  this  work  to  stationary  bifurcation.  These  results  were  successfully 
applied  to  axial  flow  compressors  [39],  [40],  to  aircraft  stall  at  high  angle-of- 
attack  [8]  and  to  voltage  collapse  [54],  In  [54]  and  [3],  the  issue  of  relocating  a 
bifurcation  while  at  the  same  time  stabilizing  the  bifurcated  solution  has  been 
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addressed.  In  [9],  the  control  design  result  of  [1]  was  extended  in  a  new  dynamic 
feedback  structure  incorporating  a  washout  filter.  The  advantage  gained  by 
employing  washout  filters  in  the  feedback  loop  is  the  preservation  of  all  system 
equilibria. 

Mehra  [45]  and  Melira,  Kessel,  and  Carroll  [46]  investigated  the  third  category 
of  bifurcation  control,  namely  global  removal  of  bifurcations  by  state  feedback. 
Their  results  can  apply  only  to  stationary  bifurcations,  since  it  is  obtained  by 
appealing  to  a  global  implicit  theorem.  These  publications  were  the  first  on  any 
aspect  of  bifurcation  control. 

Even  with  significant  control  authority,  it  might  not  be  possible  to  remove 
bifurcations  altogether.  In  fact,  much  work  in  bifurcation  control  deals  with 
simply  delaying  the  bifurcation  as  much  as  possible  and/or  stabilizing  the  bifur¬ 
cated  solution  instead  of  eliminating  it.  Applying  control  to  delay  or  eliminate  a 
bifurcation  might  be  costly  and  unnecessary  if  the  system  normally  operates  far 
from  conditions  leading  to  bifurcation.  For  such  cases,  we  may  want  to  apply 
control  only  when  the  system  approaches  bifurcation.  In  some  cases,  an  accu¬ 
rate  model  of  the  system  is  not  available  for  control  design,  or  physical  means  of 
actuation  might  not  be  sufficient  for  controlling  a  bifurcation.  For  these  cases,  it 
is  imperative  that  the  system  (or  system  operator)  recognizes  during  operation 
that  the  system  is  approaching  a  bifurcation.  Once  such  a  condition  is  detected, 
appropriate  measures  can  be  taken. 

These  considerations  motivate  us  to  develop  monitoring  techniques  that  can 
provide  some  type  of  warning  signal  that  the  system  is  close  to  a  bifurcation. 
Wiesenfeld  [57]  has  shown  that  the  power  spectrum  of  a  measured  output  exhibits 
distinguishing  features  near  a  bifurcation  for  systems  normally  operating  along  a 
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periodic  solution.  For  many  engineering  systems,  the  normal  operating  condition 
is  an  an  equilibrium  point  rather  than  a  periodic  solution.  Thus,  we  extend  this 
work  to  systems  operating  at  an  equilibrium  point.  This  will  form  a  basis  for 
our  design  of  monitoring  systems  that  exhibit  precursors  of  bifurcation  for  a 
nonlinear  system  operating  at  an  equilibrium  point. 

It  should  be  noted  that  the  term  “precursor”  has  also  been  used  in  the  liter¬ 
ature  on  stall  of  axial  flow  compressors  [26],  [21]  and  instabilities  in  combustion 
systems  [52],  The  precursors  that  have  been  suggested  for  these  applications 
differ  markedly  in  spirit  to  what  is  sought  in  this  dissertation.  In  the  cited  refer¬ 
ences,  the  initiation  of  an  instability  is  detected  by  examining  the  time  signal  as 
it  begins  to  depart  from  the  nominal  operating  state.  In  contrast,  in  this  work 
we  attempt  to  obtain  warning  signals  for  instabilities  that  provide  an  indication 
of  nearness  to  instability  as  opposed  to  an  indication  of  the  onset  of  instability. 
The  implications  for  system  operability  of  having  such  precursors  available  are 
significant. 

The  remainder  of  the  dissertation  proceeds  as  follows.  Some  basic  theoretical 
results  of  the  bifurcations  and  related  mathematical  tools  will  be  discussed  in 
Chapter  2.  The  basic  local  codimension  one  bifurcations  are  presented. 

In  Chapter  3,  the  effects  of  noise  on  systems  exhibiting  bifurcation  are  inves¬ 
tigated.  First,  we  review  the  effect  of  noise  on  systems  operating  along  a  periodic 
solution  and  recall  why  the  power  spectrum  can  be  used  as  a  useful  precursor 
for  bifurcation  of  a  periodic  solution.  Then,  we  extend  this  result  to  bifurcations 
from  an  equilibrium  point. 

The  precursor  of  Chapter  3  are  very  useful  for  detecting  Hopf  bifurcation, 
but  much  less  so  for  detecting  stationary  bifurcation.  Chapter  4  addresses  how 
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to  make  stationary  bifurcation  more  detectable  by  a  power  spectrum  precursor. 
A  new  augmented  state  is  presented  as  a  solution.  Theorem  that  states  a  new 
augmented  state  transformation  renders  stationary  bifurcation  into  Hopf  bifur¬ 
cation  without  much  knowledge  on  system  dynamics  is  shown.  The  effects  of  a 
augmented  states  transformation  on  a  nonlinear  system  are  considered  in  this 
chapter.  Stability  of  bifurcated  periodic  solution  due  to  a  transformation,  effect 
of  transformation  on  a  system  experiencing  Hopf  bifurcation,  and  reduced  order 
transformation  are  also  considered.  Extension  of  an  augmented  state  transfor¬ 
mation  to  singularly  perturbed  systems  is  investigated.  Finally,  we  suggest  a 
another  augmented  state  transformation  to  relax  some  of  assumptions  which  we 
have  made  in  a  previous  transformation. 

In  Chapter  5,  we  use  the  axial  compression  system  as  a  example.  Bifurca¬ 
tion  analysis  of  an  axial  compression  system  model  is  performed.  The  model 
undergoes  both  stationary  and  Hopf  bifurcation  as  a  parameter  varies.  The 
compression  model  provides  examination  of  a  suggested  precursor  in  a  practical 
engineering  system.  We  perform  numerical  simulation  on  axial  compressor  under 
mild  assumption  on  the  modeling  part. 

Conclusions  and  future  research  directions  are  given  in  Chapter  6. 
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Chapter  2 


Preliminaries 


In  this  chapter,  we  review  some  basic  concepts  and  results  from  bifurcation 
theory,  mulitilinear  functions,  and  random  process  theory  which  will  be  used  in 
subsequent  chapters.  The  discussion  follows  [2],  [53]  and  [28]. 

2.1  Local  Bifurcations 

In  this  section,  we  recall  the  statement  of  well  known  local  bifurcation  theo¬ 
rems.  The  term  bifurcation  refers  to  qualitative  changes  in  the  phase  portraits 
of  dynamical  systems  occurring  with  slight  variation  in  the  system  parameters. 
There  are  many  types  of  bifurcation.  Here,  we  have  particular  interest  in  local 
bifurcations,  i.e.  bifurcations  in  the  neighborhood  of  an  equilibrium  point.  The 
nominal  operating  condition  of  an  engineering  system  can  often  be  taken  to  be 
an  equilibrium  point. 

2.1.1  Low-dimensional  example 

Originally,  Poincare  used  the  term  bifurcation  to  describe  the  splitting  of  equilib¬ 
rium  solutions  for  a  family  of  differential  equations.  Bifurcations  involving  only 


5 


equilibrium  points  are  known  as  stationary  or  static  bifurcations.  There  are  also 
bifurcations,  such  as  Hopf  bifurcation,  which  involve  both  an  equilibrium  point 
and  a  periodic  solution. 

Consider  a  system 

x  =  f(x,n)  (2.1) 

where  x  G  Rn  is  the  system  state  and  //  G  Rk  denotes  a  /.'-dimensional  parameter; 
k  can  be  any  positive  integer.  In  this  work,  we  limit  k  to  be  1  so  that  fi  is  a 
scalar.  The  equilibrium  solutions  are  given  by  the  solutions  of  the  equation 
f(x ,  //)  =  0.  By  the  implicit  function  theorem,  as  n  varies,  these  equilibria 
are  smooth  function  of  //  as  long  as  Dxf ,  the  Jacobian  matrix  of  /(.r,  n)  with 
respect  to  x  evaluated  at  the  equilibrium,  does  not  have  a  zero  eigenvalue.  An 
equilibrium  point  for  a  given  parameter  value  is  called  a  “stationary  bifurcation 
point”  if  two  or  more  equilibria  join  at  that  point.  A  necessary  condition  for 
(xo,Ho)  to  be  a  stationary  bifurcation  point  is  that  the  Jacobian  Dxf  has  a 
zero  eigenvalue  for  x  =  x0,  //  =  n o-  Moreover,  it  is  also  true  that  a  necessary 
condition  for  a  local  bifurcation  any  kind  to  occur  is  that  Dxf  have  at  least 
one  eigenvalue  with  a  zero  real  part.  For  example,  a  pair  of  purely  imaginary 
eigenvalues  typically  results  in  a  Hopf  bifurcation. 

Bifurcations  are  often  classified  according  to  their  codimension.  The  codi¬ 
mension  number  is  the  minimum  number  of  parameters  that  are  needed  to  make 
the  description  of  the  dynamics  possible,  in  the  vicinity  of  the  singularity  [47] . 
A  rigorous  definition  of  codimension  can  be  found  in  [27]  and  [38].  After  an 
appropriate  linear  coordinate  transformation,  Dx  f  can  be  represented  in  block- 
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diagonal  form 


DJ  = 


f  Ac  0  X 


(2.2) 


V 


J 


0  Aa 

where  Ac  is  the  Jordan  block  corresponding  to  the  critical  eigenvalues  (i.e.  those 
with  zero  real  part)  and  As  involves  the  remaining  eigenvalues.  Bifurcations 
from  an  equilibrium  of  codimension  one  and  two  fall  into  one  of  the  following 
situations. 

Codimension  1  bifurcation: 


1.  Ac  =  0,  a  scalar. 


2.  Ac  is  2  x  2  and  has  a  pair  of  pure  imaginary  eigenvalues. 


Codimension  2  bifurcation: 


1.  Ac  is  2  x  2  and  is  nondiagonalizable  with  a  double  zero  eigenvalue,  that  is 

toil 


1, 


v 0  V 


(2,3) 


2.  Ac  is  2  x  2  and  is  diagonalizable  with  a  double  zero  eigenvalue,  that  is 

I  0  0  1 


1,  = 


V°  °/ 


(2.4) 


3.  Ac  is  3  x  3  and  has  one  zero  eigenvalue  and  one  pair  of  pure  imaginary 
eigenvalue,  that  is 


Ac  = 


^  0  -ay  0  ^ 


ay  0  0 


0  0  0 


(2,5) 
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4.  Ac  is  4  x  4  and  has  two  pairs  of  pure  imaginary  eigenvalues,  that  is 


/ 


0 


—co  i  0 


Ac  = 


u)  i  0  0  0 

0  0  0  — uj‘2 


y  0  0  ui‘2  0  J 


(2.6) 


By  employing  what  is  known  as  center  manifold  reduction  [38]  and  a  normal 
form  transformation,  system  (2.1)  can  be  reduced  to  a  lower  order  simplified 
system  called  the  normal  forms.  The  normal  form  preserves  the  qualitative 
properties  of  the  solution  near  bifurcation.  Analyzing  the  dynamics  of  normal 
forms  yields  a  qualitative  picture  of  the  solution  for  each  type  of  bifurcation. 
The  normal  forms  of  codimension  one  bifurcations  are  summarized  as  follows: 
i)  Saddle-node  bifurcation:  The  normal  form  is  given  by 


X  =  —  X2 


(2.7) 


where  x  is  a  scalar  variable.  Equilibrium  solutions  exist  only  for  [i  >  0  and  are 
given  by  x  =  ±yfji.  The  branch,  x  =  y//!,  is  stable  while  the  other  branch, 
x  =  —y/Ji,  is  unstable. 

ii)  Transcritical  bifurcation:  The  normal  form  is  given  by 

x  =  fix  —  x2  (2.8) 


where  x  is  a  scalar  variable.  The  nominal  equilibrium  point  is  the  origin  for  all 
values  of  //,.  The  bifurcated  equilibrium  solution  is  x  =  //.  and  exists  for  both 
//  >  0  and  //  <  0.  For  //  >  0  (resp.  //  <  0),  the  bifurcated  branch  is  stable  (resp. 
unstable). 

iii)  Pitchfork  bifurcation:  The  normal  form  (for  the  supercritical  case)  is 
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given  by 


x  =  px  ±  x3  (2.9) 

where  x  is  a  scalar  variable.  Again,  the  origin  exists  as  the  nominal  equilibrium 
point  for  all  value  of  //,.  There  are  two  bifurcated  equilibrium  branches.  For  “+” 
case  in  (2.9),  these  branches  are  x  =  for  //,  >  0,  and  they  are  both  stable. 
The  opposite  is  true  in  the  case. 

iv)  Hopf  bifurcation:  The  normal  form  is  two  dimensional  and  given  by 

x  =  —  y  +  x([i  ±  (x2  +  y2)) 

y  =  -X  +  y(ii±(x2  +y2))  (2.10) 

where  x  and  y  are  scalar  variables.  The  associated  bifurcated  solution  is  a  non¬ 
constant  periodic  trajectory.  In  the  “-’’case  in  (2.10),  the  periodic  solution  occurs 
for  //  >  0  and  is  stable,  while  in  the  “+”  case  it  exists  for  //,  <  0  and  is  unstable. 


2.1.2  Bifurcation  theorems 

We  have  the  following  codimension  one  bifurcation  theorems  for  the  system  (2.1). 

Theorem  2.1  (Stationary  Bifurcation  Theorem)  Suppose  f  (x,  p)  of  system  (2.1) 
is  sufficiently  smooth  with  respect  to  both  x  and  p,  and  /( 0,  p)  =  0  for  all  p,  and 
the  Jacobian  of  f(x,p),  A ^  :=  Dxf(0,p),  possesses  a  simple  eigenvalue  X(p) 
such  that  at  the  critical  parameter  value  pc  =  0, 

A'(0)==^|(_^0  (2.11, 

and  all  the  remaining  eigenvalues  of  A0  have  strictly  negative  real  part.  Then: 
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i)  there  is  an  eo  >  0  and  a  function 


p{e)  —  p\e  +  ^e2  +  0(e 3)  (2.12) 

such  that  if  p i  ^  0;  there  is  a  nontrivial  equilibrium  x(p)  near  x  =  0  for 
each  e  G  {[— eo,  0)  U(0,  eo)};  if  hi  =  0  and  p2  >  0  (resp.  p  <  0),  there 
are  two  equilibrium  points  x±(p)  near  x  =  0  for  each  p  G  (0,  eo]  (resp. 
p  G  [— eo,  0)J 

ii)  Exactly  one  eigenvalue  [3(e)  of  the  Jacobian  evaluated  with  respect  to  each 
of  the  nontrivial  equilibrium  points  in  (i)  approaches  0  as  e  J,  0  and  it  is 
given  by  a  real  function 

/ 3(e)  =  (3\e  +  +  0(0^)  (2.13) 

The  coefficient  0\  of  this  function  satisfies  j3\  =  —  A'(0)^i.  The  nontrivial 
equilibrium  x_  (resp.  x+)  is  stable  (resp.  unstable)  if  f3\e  <  0  and  is 
unstable  (resp.  stable)  if  (3\t  >  0.  Nevertheless,  the  bifurcation  point  itself 
is  unstable.  If  (3 1  =  0,  then  02  =  — 2A'(0 )p2,  and  the  nontrivial  equilibria 
are  asymptotically  stable  if  f32  <  0  but  are  unstable  if  (32  >  0. 

The  assumptions  of  the  theorem  above  are  not  generic  for  the  case  of  a  single 
parameter  p.  A  result  that  is  generic  Theorem  2.2  below,  which  gives  conditions 
for  saddle  node  bifurcation.  In  this  bifurcation,  nominal,  stable  equilibrium 
merges  with  another,  unstable  equilibrium. 

To  state  the  theorem  on  saddle  node  bifurcation,  we  consider  system  (2.1) 
where  /  is  sufficiently  smooth  and  /( 0,  0)  =  0.  Express  the  expansion  of  f(x,  p) 
in  a  Taylor  series  about  x  =  0,  p  =  0  in  the  form 

f(x,  p)  =  Ax  +  bp  +  Q(x,  x)  +  •  •  •  (2-14) 
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Note  that  A  =  Dxf( 0,  0)  is  simply  the  Jacobian  matrix  of  /  at  the  origin  for 

H  =  0. 


(SN1)  The  Jacobian  A  possesses  a  simple  zero  eigenvalue. 

(SN2)  lb  ^  0  and  lQ(r,r )  ^  0  where  r  (resp.  1)  is  the  right  column  (resp.  row) 
eigenvector  of  the  Jacobian  A  associated  with  the  zero  eigenvalue,  with  r 
and  l  normalized  by  setting  the  first  component  of  r  to  1  and  then  choosing 
l  such  that  Ir  =  1. 

Theorem  2.2  (Saddle  Node  Bifurcation  Theorem)  If  (SN1)  and  (SN2)  hold, 
then  there  is  an  eo  >  0  and  a  function 

/t(e)  =  /^g2  +  0(c3)  (2.15) 

such  that  //2  ^  0  and  for  each  e  G  (0,e0]>  (%-l)  has  a  nontrivial  equilibrium  x(e) 
near  0  for  //  =  //(e).  The  equilibrium  point  x  =  0  is  unstable  at  bifurcation,  i.e., 
for  ji  =  0. 

Suppose  that  instability  of  the  nominal  equilibrium  x(p0)  is  the  result  of  a 
pair  of  eigenvalues  of  the  system  linearization  crossing  the  imaginary  axis  in  the 
complex  plane.  Then,  as  is  well  known  [34],  [42],  [47],  generically  it  will  be  the 
case  that  a  small  amplitude  periodic  orbit  of  (2.1)  emerges  from  the  equilibrium 
:7’o  (//)  at  the  critical  parameter  value.  We  have  the  following  theorem  for  the 
case. 

Theorem  2.3  (Cr -Hop f  Bifurcation  Theorem)  [34]  Suppose  the  system  (2.1) 
satisfies  the  following  conditions: 

i)  fn( 0)  =  0  for  //  in  an  open  interval  containing  0,  and  0  G  Rn  is  an  isolated 
equilibrium  point  of  f . 
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ii)  All  partial  derivatives  of  the  components  flp  of  the  vector  f  of  orders  l  <  r 
r  >  4  (including  the  partial  derivative  with  respect  to  p)  exist  and  are 
continuous  in  x  and  p,  in  a  neighborhood  of  (0,  0)  in  Rn  x  R 1  space. 

Hi)  :=  Dxf( 0,  p)  has  a  complex  conjugate  pair  of  eigenvalues  X  and  A  such 
that  X(p)  =  a(p)  +  jca(p),  where  w o  :=  w( 0)  >  0.  «(0)  =  0,  and 

dev 

«,(°):=— |/J=°^0  (2.16) 

iv)  The  remaining  eigenvalues  of  A0  have  strictly  negative  real  parts. 

Then: 

i)  There  exist  an  ep  >  0  and  Cr~l  function 

p(t)  =  X]  +  0(e?  x)  (0  <  e  <  ep)  (2-17) 

2=1 

such  that  for  each  e  G  (0,ep)  there  exists  a  nonconstant  periodic  solution 
p€(t)  with  period 

27T  ^ 

X(e)  =  — [1  +  X]  r2 if11}  +  0(A  x)  (0  <  e  <  ep)  (2.18) 

wo  i=i 

occurring  for  p  =  p(e). 

ii)  There  exists  a  neighborhood  r)  of  x  =  0  and  an  open  interval  A  containing 
0  such  that  for  any  p  G  A,  the  only  nonconstant  periodic  solutions  that  lie 
in  7)  are  members  of  the  family  pe(t) . 

in)  Exactly  two  of  the  Floquet  exponents  of  p£(t)  approach  0  as  e  f  0.  One  is 

0  identically,  and  the  other  is  a  C”'-1  function 

[1=2] 

/3(e)  =  X]  /^e2*  +  ^(e?  X)  (0  <  e  <  ep)  (2.19) 

2=1 
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The  periodic  solution  p€(t)  is  orbitally  asymptotically  stable  if  (3 (e)  <  0, 
and  is  unstable  if  /3(e)  >  0.  If  there  is  a  first  nonvanishing  coefficient  p2k, 
then  the  first  nonvanishing  coefficient  in  (2.19)  is  given  by 

I32k  =  -2\'(0)p2k  (2.20) 

Moreover,  there  is  then  an  e\  G  (0,ep)  such  that 

sgn[/3( p)]  =  sgn[/32]  (2.21) 

for  p  G  {^|0  <  p/p2 k  <  p(ef)/pik)-  Here,  sgn  denotes  the  sign  of  a  real 
number. 

2.1.3  Singularly  perturbed  Hopf  bifurcation 

Since  we  are  also  going  to  consider  the  case  of  Hopf  bifurcation  under  singular 
perturbation  in  this  thesis,  discussion  of  a  basic  theorem  on  this  topic  is  in  order. 
The  next  theorem  gives  conditions  that  ensure  persistence  of  Hopf  bifurcation 
under  singular  perturbation  of  a  dynamical  system.  Let  the  full  system  be 

x  =  f(x.  y,  p,  e) 

ey  =  g(x,y,p,e)  (2.22) 

where  x  G  Rn ,  y  G  Rm ,  //.  e  G  R  and  e  is  small  but  nonzero.  Here,  p  is  the  bifur¬ 
cation  parameter  and  a  e  is  the  singular  perturbation  parameter.  The  associated 
reduced  system  ,  obtained  by  formally  setting  e  =  0  in  (2.22),  is 

x  =  f(x,  y,p,  0) 

0  =  g(x,  y,  p,  0)  (2.23) 
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Let  J(fjL,e)  denote  the  Jacobian  matrix 


(  DJ  Dyf  ^ 

v  e~1Dxg  e~lDyg  y 


(2.24) 


of  the  full  system  evaluated  at  the  equilibrium  point  Xq (/jl,  e),  yo(g,  e)  of  interest. 
Introduce  the  following  assumptions 


(51)  f.9  are  Cr  (r  >  5)  in  x,  y,  ji.  e  in  a  neighborhood  of  (x0,  y0, 0,  0). 

(52)  det  Dyg(x0,  y0,  0,  0)  ^  0. 

(53)  The  reduced  system  (2.23)  undergoes  a  Hopf  bifurcation  from  the  equilib¬ 
rium  m0  =  (%o,y o)  for  the  critical  parameter  value  //  =  0. 

(54)  No  eigenvalue  of  D2g{x 0,  yo,  0,  0)  has  zero  real  part. 

Under  conditions  (S1)-(S4),  Abed  [1]  has  proved  the  following  theorem. 

Theorem  2.4  (Persistence  Under  Singular  Perturbation)  Let  (S1)-(S4)  above 
hold.  Then,  there  is  an  e0  >  0  and  for  each  e  €  (0,  e0]  the  full  system  (2.22) 
undergoes  a  Hopf  bifurcation  at  an  equilibrium  near  m0  for  a  critical  pa¬ 

rameter  value  pLec  near  0.  We  also  have 

lim  —  Re(\i(fiec,  e))  =  et^O)  (2.25) 

where  e),  \i(n,  e)  are  the  complex- conjugate  eigenvalues  of  J(n,e)  which 
cross  the  imaginary  axis  for  g  =  gec,  and  where  a(g)  is  the  real  part  of 
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2.2  Multilinear  Functions  and  the 


Fredholm  Alternative 

2.2.1  Multilinear  functions 

Multivariable  Taylor  series  expansions  are  used  in  the  derivation  of  local  stability 
conditions  for  bifurcations.  The  notation  of  multilinear  functions  is  convenient  in 
dealing  with  multivariable  Taylor  series.  The  following  are  some  basic  definitions 
and  properties  of  mulitilinear  functions. 

Definition  2.1  Let  V\ .  V> ■ . . . ,  14  and  W  be  vector  spaces  over  the  same  field. 
A  map  fi  :  \'\  x  V2  x  •  •  •  x  14  e- ?•  W  is  said  to  be  multilinear  (or  k-linear)  if  it  is 
linear  in  each  of  its  variables.  That  is  ([15],  p.  76),  for  arbitrary  v\  if  E  \ 4  i  = 
1 , ...  ,k,  and  for  arbitrary  scalars  a ,  a,  we  have 

fiiv1,  •  •  • ,  avl  +  civ\  ■  ■  ■ ,  vk)  =  affv1,  •  •  •  ,v1,  -  ■  ■ ,  vk) 

+  dfi{v\---,v\---,vk)  (2.26) 

We  refer  to  k  as  the  degree  of  the  multilinear  function  fi.  In  particular,  multi¬ 
linear  functions  of  degree  two,  three,  and  four  are  referred  as  bilinear,  trilinear, 
and  tetralinear  functions ,  respectively. 

In  the  sequel,  we  shall  deal  exclusively  with  mulitilinear  functions  whose 
domain  is  the  product  space  of  k  identical  vector  spaces  V\  =  Vi  =  •  •  •  =  14  =  V. 
For  such  multilinear  functions,  we  have  the  following  notion  of  symmetry. 

Definition  2.2  A  k-linear  function  fi  :  V  x  V  x  •  •  •  x  V  — >■  W  is  symmetric  if, 
for  any  if  E  V,  i  =  1, . . . ,  k  the  vector 

f’(v\v2,---,vk)  (2.27) 
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is  invariant  wider  arbitrary  permutations  of  the  argument  vectors  vl . 

With  an  arbitrary  multilinear  function  'W.  we  associate  a  symmetric  multi¬ 
linear  function  yy  by  the  symmetrization  operation  ([15],  pp.  88-89).  Given 
a  multilinear  function  ffx1 ,  x2 ,  ■  ■  ■ ,  xk),  define  a  new  (symmetric)  multilinear 
function  yy  as  follows: 

•  •  •  ,xk)  :=  y  (2.28) 

where  the  sum  is  taken  over  the  k\  permutation  of  the  integer  1,  2, . . . ,  k. 

An  important  property  of  homogeneous  functions  represented  in  terms  of 
multilinear  fictions  is  given  next. 

Proposition  2.1  Let  rtp  :  ( Rm)k  — >■  Rm  be  a  symmetric  k-linear  function.  For 
any  vector  v  G  Rn, 

DiKv,  G  v  =  klH'>h  nr--,  V-  v)  (2.29) 

2.2.2  The  Fredholm  Alternative 

Next,  we  consider  the  Fredholm  Alternative.  Consider  the  system  of  linear  equa¬ 
tions 

Ax  =  b  (2.30) 

where  A  is  a  real  n  x  n  matrix  and  b  G  Rn.  If  the  coefficient  matrix  A  is 
nonsingular,  then  (2.30)  has  a  unique  solution,  given  by  A~lb.  Existence  of  solu¬ 
tions  for  cases  in  which  A  is  singular  is  the  subject  of  the  Fredholm  Alternative. 
Below,  we  first  present  this  result  for  the  general  case  of  a  singular  coefficient 
matrix  A,  and  then  employ  it  for  the  particular  case  in  which  A  has  a  simple 
zero  eigenvalue. 
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Theorem  2.5  Let  N(A )  be  k- dimensional,  with  basis  r1,  ■  ■  ■ ,  rk ,  and  dual  basis 
l\  ■  ■  ■ ,  lk.  Then  (2.30)  has  at  least  one  solution  in  R"  if  and  only  if  llb  =  0  for 
i  =  1  ,...,k.  Moreover,  in  such  a  case,  the  general  solution  of  (2.30)  has  the 
representation 

k 

x  =  x°  +  ^2  (2.31) 

2=1 

where  a:0  is  any  particular  solution  of  (2.30)  and  the  are  arbitrary  real  scalar. 

Suppose  now  that  A  has  a  simple  zero  eigenvalue.  Let  r  and  l  denote  right 
and  left  eigenvectors  of  A,  respectively,  corresponding  to  the  zero  eigenvalue, 
and  require  that  these  be  chosen  to  satisfy  Ir  =  1.  Under  these  conditions,  the 
Fredholm  Alternative  asserts  that  (2.30)  has  a  solution  if  and  only  if  lb  =  0. 
Moreover,  the  Fredholm  Alternative  also  implies  that,  if  (2.30)  has  a  solution  x° , 
then  the  totality  of  solutions  is  given  by  the  one-parameter  family  x  =  .r0  +  ar 
where  a  G  R  is  arbitrary.  The  solution  is  rendered  unique  upon  imposing  a 
normalization  condition  which  specifies  the  value  of  lx. 

Introduce  subspaces  Ec ,  Es  G  Rn  as  follows:  Ec  is  the  one-dimensional  sub¬ 
space 

Ec  :=  span{r}  (2.32) 

and  Es  is  the  (n  —  l)-dimensional  subspace 

Es  :=  {x  e  Rn\lx  =  0}  (2.33) 

From  the  foregoing,  we  have  in  particular  that  if  lb  =  0  then  system  Ax  = 
b ,  lx  =  0  has  a  unique  solution.  Equivalently,  (2.30)  has  a  unique  solution  in  Es 
for  any  vector  b  G  Es .  This  proves  that  the  restriction  A\Es  of  the  linear  map 
defines  an  invertible  (one-to-one  and  onto)  map.  In  the  next  result,  we  exhibit 
the  unique  solution  which  lies  in  Es  of  the  system  Ax  =  0,  lx  =  0. 
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Proposition  2.2  Suppose  A  has  a  simple  zero  eigenvalue.  Then  the  unique 
solution  of  Ax  =  b ,  lx  =  0  given  that  lb  =  0  is 

x  =  (ATA  +  lTl)~1ATb  (2.34) 

Therefore,  the  inverse  of  the  restricted  map  A\Es  exists  and  is  given  by 

(AIes)-1  =  ( AtA  +  lTl)~1AT  (2.35) 

2.3  Stability  of  Bifurcated  Solutions 

Abed  and  Fu  [5], [6]  derived  formulas  for  determining  the  local  stability  of  bifur¬ 
cated  solutions.  We  review  their  method  here  since  we  will  employ  it  later  to 
study  the  stability  of  the  bifurcated  solution  of  a  system  system  monitored  to 
detect  possible  nearness  to  bifurcation. 

Consider  a  one-parameter  family  of  nonlinear  autonomous  systems 

x  =  f{x,qi)  (2.36) 

where  x  G  R"  is  the  vector  state  and  //  is  real-valued  parameter.  Let  f{x,g)  be 
sufficiently  smooth  in  x  and  qi  and  let  x Q  be  the  nominal  equilibrium  point  of 
the  system  as  a  function  of  the  parameter  //,. 

2.3.1  Stability  calculation  for  stationary  bifurcation 

First,  we  consider  the  case  of  stationary  bifurcation.  The  next  hypothesis  is  in 
force  throughout  this  section. 

(S)  The  Jacobian  matrix  of  system  (2.36)  at  the  equilibrium  x®  has  a  simple 
zero  eigenvalue  Ai (n)  with  A^O)  ^  0,  and  the  remaining  eigenvalues  lie  in 
the  open  left  half  of  the  complex  plane  for  jic  =  0. 
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Theorem  2.1  asserts  that  hypothesis  (S)  implies  a  stationary  bifurcation  from 
Xq  at  /i  =  0  for  (2.36).  That  is,  a  new  equilibrium  point  bifurcates  from  x®  at 
//  =  0.  Recall  that  near  the  point  (a^O)  of  the  (n+  l)-dimensional  (x,  fi) -space, 
there  exists  a  parameter  e  and  a  locally  unique  curve  of  critical  points  (rc(e),  p(e)), 
distinct  from  x(]jA  and  passing  through  (.Xq,  0),  such  that  for  all  sufficient  small  |e|, 
x(e)  is  an  equilibrium  point  of  (2.36)  when  //,  =  p(e). 

The  parameter  e  may  be  chosen  such  that  x(e),n(e)  are  smooth.  The  series 
expansions  of  x(e),n(e)  can  be  written  as 

P'(e)  =  fiie  +  P2^2  +  ■  ■  ■  (2.37) 

x(e)  =  x0^  +  x\e X2t2  +  ■  ■  ■  (2.38) 

If  Hi  7^  0,  the  system  undergoes  a  transcritical  bifurcation  from  x °  at  //  =  0. 
That  is,  there  is  a  second  equilibrium  point  besides  x(]jf  for  both  positive  and 
negative  values  of  //  with  |p|  small.  If  ji\  =  0  and  //•;  ^  0,  the  system  undergoes 
a  pitchfork  bifurcation  for  \n\  sufficiently  small.  That  is,  there  are  two  new 
equilibrium  points  existing  simultaneously,  either  for  positive  or  for  negative 
values  of  //  with  |//|  small.  The  new  equilibrium  points  have  an  eigenvalue  (3(e) 
which  vanishes  at  //  =  0.  The  series  expansion  / 3(e)  is  given  by 

/ 3(e)  =  p\e  +  p2e2  +  •  •  •  (2.39) 

with 

Pi  =  -^iA'(O)  (2.40) 

and,  in  case  /?i  =  0.  (3>  is  given  by 

P2  =  -2^(0)  (2.41) 
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The  stability  coefficients  (3 1  and  p2  can  be  determined  solely  by  eigenvector 
computations  and  the  coefficients  of  the  series  expansion  of  the  vector  field. 
System  (2.36)  can  be  written  in  the  series  form 

x  =  L^x  +  Qn(x,  x)  +  Cfl{x,  x,  £)  +  ••• 

=  LqX  pLix  T  f  /. t  T  •  •  • 
p-Qo(x,x)  +  nQi{x,x)  4 - 

+Cq(x,  x.  x)  I  •  •  •  (2.42) 

where  x  =  x  —  x®,  LfJ,  Li,  L2  are  n  x  n  matrices,  Q^ix.  x).  Qo(x-  x).  Q\  (x.  x) 
are  vector-valued  quadratic  forms  generated  by  symmetric  bilinear  forms,  and 
C^ix,  x ,  x),  Co(x,  x,  x)  are  vector-valued  cubic  forms  generated  by  symmetric  tri- 
linear  forms. 

By  assumption,  the  Jacobian  matrix  L0  has  only  one  simple  zero  eigenvalue 
with  the  remaining  eigenvalues  stable.  Denote  by  /  and  r  the  left  (row)  and 
right  (column)  eigenvectors  of  the  matrix  L0  associated  with  the  simple  zero 
eigenvalue,  respectively,  where  first  component  of  r  is  set  to  be  1  and  the  left 
eigenvector  l  is  chosen  such  that  Ir  =  1.  It  is  well  known  that 

A'(0)  =  lLxr  (2.43) 

A  stability  criterion  for  the  bifurcated  equilibria  of  system  (2.42)  is  given  in  the 
following  two  lemmas. 

Lemma  2.1  The  two  bifurcated  equilibrium  points  of  (2.f2)  near  x Q  for  qi  near 
0,  which  appear  only  for  p  >  0  (resp.  p  <  0)  when  ILp'  >  0  (resp.  IL^r  <  0), 
are  asymptotically  stable  if  pi  =  0  and  [3  <  0,  and  unstable  if  pi  =  0  and  (3  >  0. 
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Here, 


(2.44) 

(2.45) 


Pi  =  lQo(r,r) 

P2  =  2/ { 2Q0  ( r,  x2 )  + C0  ( r,  r,  r ) } 
with  x2  satisfies  the  following  equation: 

L0x2  =  -Q0{r,  r)  (2.46) 


Lemma  2.2  Suppose  the  value  of  pi  given  in  (2-44)  negative.  Then  the  bi¬ 
furcated  solution  occurring  for  for  fi>  0  (resp.  fi  <  0)  is  asymptotically  stable 
when  ILyr  >  0  (resp.  lL\r  <  0 ). 

The  criterion  given  in  Lemma  2.1  corresponds  to  pitchfork  bifurcation,  with 
the  one  in  Lemma  2.2  is  for  transcritical  bifurcation. 

2.3.2  Stability  calculation  for  Hopf  bifurcation 

Now  consider  system  (2.36)  under  for  the  following  hypothesis,  which  implies 
occurrence  of  Hopf  bifurcation 

(H)  The  Jacobian  matrix  of  system  (2.36)  at  the  equilibrium  x °  has  a  pair  of 
pure  imaginary  eigenvalues  Xi(fic)  =  icoc  and  Xi(fic)  =  —hvc  with  u>c  ^  0, 
the  transversality  condition  dReWvc)\  ^  o  is  satisfied,  and  all  the  remaining 
eigenvalues  lie  in  the  open  left  half  complex  plane. 

Theorem  2.3,  the  Hopf  bifurcation  theorem  asserts  the  existence  of  a  one- 
parameter  family  p(,  0  <  e  <  eo  of  nonconstant  periodic  solutions  of  system 
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(2.36)  emerging  from  x  =  x®c  at  the  parameter  value  jic  for  sufficiently  small 
\p  —  nc\.  Exactly  one  of  the  characteristic  exponents  of  pe  governs  the  asymptotic- 
stability  and  is  given  by  a  real,  smooth  and  even  function 

/3(e)  =  /32e2  +  /3 4c4  +  •  •  •  (2.47) 

Specifically,  pe  is  orbitally  stable  if  /3(e)  <  0  but  is  unstable  if  /3(e)  >  0.  Gener¬ 
ic-ally  the  local  stability  of  the  bifurcated  periodic  solution  p,  is  decided  by  the 
sign  of  the  coefficient  /32.  Note  the  sign  of  /32  also  determines  the  stability  of 
the  critical  equilibrium  point  x°^c.  An  algorithm  for  computing  the  stability 
coefficient  /32  is  given  as  follows: 

Stepl  Express  (2.36)  in  the  Taylor  series  form  (2.42).  Let  r  be  the  right  eigen¬ 
vector  of  L0  corresponding  to  eigenvalue  icoc  with  the  first  component  of 
r  equal  to  1.  Let  l  be  the  left  eigenvector  of  L0  corresponding  to  the 
eigenvalue  iuc.  normalized  such  that  Ir  =  1. 

Step2  Solve  the  equations 

L0ci  =  ~^Qo{r,r)  (2.48) 

(2 iujcI  -  L0)b  =  ^Qo{r,r)  (2.49) 

for  a  and  b. 

Step  3  The  stability  coefficient  ,/32  is  given  by 

,/32  =  2Re[2lQ0(r,  a)  +  lQ0{r ,  b )  +  ^/C(r,  r,  r)]  (2.50) 
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2.4  Washout  Filters  in  Nonlinear  Control 


Washout  filters  are  used  commonly  in  control  systems  for  power  systems  and 
aircraft.  The  main  advantage  of  using  these  filters  is  the  resulting  robustness 
of  the  system  operating  point  to  model  uncertainty  and  to  other  control  actions 
which  may  be  used.  In  this  section,  we  give  a  brief  discussion  of  these  filters,  their 
use  in  control  of  parameterized  systems  and  especially  system  with  bifurcation. 

A  washout  filter  (or  washout  circuit)  is  a  stable  high  pass  filter  with  transfer 
function  [25] 


G(s) 


y{s)  =  s 

x(s)  s  +  a 


s  +  a 

Here,  a  >  0  is  the  inverse  of  the  filter  time  constant.  Letting 


(2.51) 


z(s)  :=  - x(s)  (2.52) 

s  +  a 

the  dynamics  of  the  filter  can  be  written  as 


z  =  x  —  az 


(2.53) 


and 

y  =  x  —  az  (2.54) 

Feedback  through  washout  filters  results  in  equilibrium  preservation  in  the 
presence  of  system  uncertainties  and  other  control  mechanisms  that  they  inher¬ 
ently  offer.  Indeed,  the  most  striking  property  of  a  system  controlled  through  a 
washout  filter  is  that  all  the  original  equilibria  are  preserved.  Equilibrium  points 
represent,  in  some  sense,  a  system’s  capability  to  perform  in  a  certain  manner 
at  steady  state.  There  are  situations  in  which  such  a  capability  should  not  be 
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altered  by  the  introduced  control  strategy,  such  as  in  the  lateral  control  design 
for  an  aircraft. 

Washout  filters  reject  steady-state  inputs,  while  passing  transient  inputs.  At 
steady- state, 


x 

a 


the  output  y  =  0,  and  the  steady-state  input  signal  has  been  washed  out. 


Consider  a  system 


(2.55) 


X  =  f(x,  u) 


(2.56) 


with 

f{x  o,0)  =0  (2.57) 

where  u  is  the  control  input  and  x$  is  an  equilibrium  point  for  the  system  with 
zero  input.  Let  the  control  input  u  be  a  function  of  y  (denote  it  u  =  h(y )),  and 
let  h  satisfy 

h(0)  =  0  (2.58) 

From  (2.53)-(2.54),  it  is  clear  that  y  vanishes  in  steady  state.  Hence 


f(xo,h{y0 ))  =  f[x  0,0)  =  0 


(2.59) 


and  xq  remains  an  equilibrium  point  of  the  closed-loop  system.  This  shows  that 
by  incorporating  a  washout  filter  in  the  feedback,  the  equilibrium  points  of  the 
original  system  are  preserved. 


2.5  Background  on  the  Power  Spectral  Density 

In  this  section,  we  review  the  definition  of  the  power  spectral  density  of  a  random 
process  and  discuss  related  computational  aspects.  We  begin  by  recalling  basic 
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definitions  from  the  theory  of  random  processes.  The  definitions  given  next 
follow  Gray  and  Davisson  [28],  and  the  discussion  of  computation  of  the  power 
spectral  density  follows  Jenkins  and  Watts  [37]. 

A  random  process  or  stochastic  process  is  an  indexed  family  of  random  vari¬ 
ables  {Xt  :  t  €  1}  defined  on  a  common  probability  space  (Q,F,P)  where  Q 
is  sample  space,  F  is  a  sigma  field  of  a  sample  space  Q.  and  P  is  a  probability 
measure  on  a  measurable  space  (Q,  F).  The  process  said  to  be  discrete-time  if  I 
is  discrete  and  continuous-time  if  the  index  set  I  is  continuous. 

Given  a  random  process  {Xt  :  t  G  /},  the  autocorrelation  function  R  x(t,  s )  : 
t,  s  G  /  is  defined  by 

Rx(t,s)  =<  XtXs  >  (2.60) 

The  autocovariance  function  l\  \(l.  s )  :  t,  s  El  is  defined  by 

Kx{t,s)  =  cov(Xt,Xs) 

=  Rx {t,  s )  —  (<  Xt  ><  Xs  >)  (2.61) 

where  cov  denotes  covariance.  Note  that  the  autocorrelation  and  autocovariance 
functions  are  equal  if  the  process  has  zero  mean,  that  is,  if  <  Xt  >=  0  for  all  t. 

The  random  process  is  said  to  be  weakly  stationary  (or  wide  sense  stationary ) 
if 

<  Xt  >=<  Xs  >  (2.62) 

for  all  t,  s  €  /  and 

Rx{t,t  +  r)  =  Rx(r)  (2.63) 

for  all  t,  t  such  that  t,  t  +  r  e  /. 

For  a  given  weakly  stationary  random  process  {Xt}  with  autocovariance  func¬ 
tion  K x in),  the  power  spectral  density  Sx(t))  is  defined  as  the  Fourier  transform 
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of  the  autocovariance  function,  that  is, 


J2kKx(k)e  j2nv<  in  discrete  time 
Sx(v)  =  \  (2-64) 

/  Kx{r)e  j27r7?T dr ,  in  continuous  time 

A  weakly  stationary  random  process  {Xt}  is  said  to  be  white  if  its  power 
spectral  density  is  a  constant  for  all  /,  that  is 

Sx(n)  =  4  (2.65) 

for  all  /. 

Let  {A",}  be  a  weakly  stationary  random  process,  and  suppose  that  Xt  is  the 
input  signal  to  a  linear  time-invariant  system  with  transfer  function  H(r) )  which 
represents  a  frequency  response  of  a  linear  time-invariant  system  to  sinusoidal 
signal  input.  If  {X,  }  has  spectral  density  Sx(t /),  then  the  output  random  process 
{!)}  has  spectral  density 

Sy(V)  =  \H(V)\2Sx(V)  (2.66) 

To  estimate  the  power  spectrum  in  an  experimental  or  computer  simulation 
setting,  there  are  two  important  factors  to  consider  [37].  First,  the  sampling- 
rate  should  be  chosen  such  that  the  estimated  power  spectrum  is  valid  within 
the  frequency  range  of  interest  0  <  7)  <  7)q.  Second,  the  length  of  the  record 
that  is  needed  to  detect  a  peak  in  the  power  spectrum  of  width  w  is  inversely 
proportional  to  w.  This  means  that  a  longer  time  history  is  needed  for  detecting  a 
peak  of  smaller  width.  For  a  fixed  length  of  the  time  record,  there  are  techniques 
for  improving  the  resolution  [37,  pp.  239-248].  The  power  spectral  densities 
given  in  the  examples  later  in  the  dissertation  are  computed  by  simulation  using 
the  Simulink  package.  In  each  case,  the  results  were  verified  by  checking  that 
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variation  of  the  record  length  and  sampling  rate  did  not  influence  the  calculated 
spectrum. 
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Chapter  3 


Noisy  Precursors  for  Nonlinear 
System  Instabilities 


In  this  chapter,  we  first  summarize  the  work  of  Wiesenfeld  [57], [58], [60]  on  prop¬ 
erties  of  nonlinear  systems  with  noise  inputs  in  the  vicinity  of  bifurcation,  and 
then  extend  these  results  both  in  scope  and  in  concept. 

The  original  motivation  for  Wiesenfeld’s  work  was  the  determination  of  the 
precise  value  of  bias  current  in  a  nonlinear  circuit  for  which  the  onset  of  a  period 
doubling  instability  occurs.  The  power  spectrum  of  a  system  variable  was  mea¬ 
sured.  Near  a  period  doubling  bifurcation,  the  power  spectrum  would  show  a  new 
peak  at  half  the  fundamental  frequency  of  the  nominal  periodic  solution  when 
period  doubling  bifurcation  occurred.  However,  substantial  broad-band  noise 
also  centered  at  half  the  fundamental  frequency  made  precise  determination  of 
the  onset  of  instability  impossible.  This  observation  motivated  the  development 
of  theory  which  could  be  applied  to  all  codimension  one  bifurcations  of  a  periodic- 
solution. 

Wiesenfeld  showed  that  the  power  spectrum  of  a  measured  output  for  a  sys¬ 
tem  operating  at  a  periodic  steady  state  and  forced  by  small  white  Gaussian  noise 
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perturbations  exhibits  a  sharply  growing  peak  near  the  fundamental  frequency 
as  the  system  nears  a  bifurcation.  Then  the  suggested  analysis  was  successfully 
applied  in  signal  amplification  area  such  as  a  resonantly  driven  silicon  p  —  n 
junction  [36]  and  a  modulated  semiconductor  laser  [59]. 

Based  on  this  observation,  we  will  introduce  the  measured  output  power 
spectrum  as  a  precursor  that  could  predict  proximity  to  bifurcation. 

3.1  Noisy  Precursor  for  System  Operating  at 
Limit  Cycle 

In  this  section, we  will  briefly  review  work  of  Wiesenfeld  on  precursors  for  bifur¬ 
cation  for  system  operating  at  a  limit  cycle.  We  do  not  give  the  detail  derivation, 
since  a  similar  derivation  will  be  used  in  next  section.  Consider  the  system 

x  =  f{x,p)  (3.1) 

where  x  G  Rn  and  //  is  a  scalar  bifurcation  parameter.  Suppose  the  system  has 
an  asymptotically  stable  T— periodic  solution  xp  for  that  changes  with  //  (T  can 
depend  on  ji) .  Thus, 

xp{t  +  T)  =  xp{t)  (3.2) 

For  small  perturbations,  the  dynamical  system  behavior  can  be  described  by 
the  linearized  equation  of  (3.1)  about  xp.  Suppose  that  the  noise  enters  as  an 
additive  forcing  term.  The  linearized  system  with  noise  forcing  term  becomes 

x  =  Df(xp)x  +  n(t)  (3.3) 
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where  x  :=  x  —  xp  and  n(t)  €  Rn  is  a  vector  white  Gaussian  noise  with  mean 
zero  and  correlation 

('ni{t)nj(t  +  r))  =  4^8{t)  (3.4) 

where  the  4>ij  are  finite  constants  for  all  i,j.  Then,  we  can  solve  equation  (3.3) 
by  using  the  fundamental  matrix  and  Floquet  theory.  A  detailed  derivation  is 
not  presented  here  since  it  is  analogous  to  the  derivation  of  the  measured  output 
power  spectrum  for  nonlinear  systems  operating  on  equilibrium  point  which  is 
presented  in  next  section. 

Using  this  solution  and  assuming  that  the  real  parts  of  the  non  critical  Flo¬ 
quet  exponents  of  the  fundamental  matrix  are  significantly  larger  (in  magnitude) 
than  real  parts  of  the  critical  Floquet  exponents  for  a  system  near  bifurcation, 
Wiesenfeld  [57]  approximated  the  solution  of  equation  (3.3)  in  terms  of  only  the 
imaginary  axis  crossing  Floquet  exponents. 

Finally,  Wiesenfeld  [57]  studied  the  power  spectrum  of  the  approximate  so¬ 
lution,  focusing  on  changes  that  might  occur  near  a  bifurcation.  Moreover,  he 
showed  that  different  types  of  bifurcation  correspond  to  different  power  spectrum 
change. 

Figure  3.1  from  [57]  shows  the  expected  effects  in  the  power  spectrum  for  each 
type  of  codimension  one  bifurcation.  This  figure  concerns  only  bifurcation  from 
a  limit  cycle,  not  from  an  equilibrium  point.  Note  that  the  power  spectrum  for 
saddle-node  or  transcritical  bifurcation  does  not  exhibit  new  peaks,  but  rather 
has  an  increasing  bump  around  the  fundamental  frequency  (defined  ^  where 
T  is  period).  Here,  fundamental  frequency  is  assumed  to  be  1  except  pitchfork 
bifurcation  that  has  0.5  fundamental  frequency.  The  other  bifurcations  such 
as  pitchfork,  Hopf  and  period  doubling  bifurcation  show  the  increasing  power 
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spectrum  peaks  at  new  locations. 


3.2  Noisy  Precursor  for  System  Operating  at 
an  Equilibrium  Point 

In  this  section,  we  consider  the  effect  of  noise  on  systems  operating  at  an  equilib¬ 
rium  point  under  condition  that  could  give  rise  to  a  bifurcation.  In  the  foregoing 
section,  we  discussed  the  effect  of  noise  in  the  analogous  situation  of  a  system 
operating  at  a  periodic  solution.  It  was  suggested  that  the  power  spectrum  could 
be  used  to  indicate  the  closeness  to  bifurcation.  Moreover,  we  saw  that  the  power 
spectrum  could  also  be  used  to  discern  the  bifurcation  type  such  as  saddle-node, 
pitchfork,  and  Hopf  bifurcation.  The  conclusions  motivate  our  work  this  section, 
where  we  study  be  the  effect  of  noise  on  the  power  spectrum  of  system  outputs 
when  an  operating  condition  is  near  bifurcation.  We  follow  closely  the  method 
used  in  [57]. 

Consider  the  autonomous  dynamic  system 

x  =  f(x,fi)  (3.5) 

where  x  E  R"  and  //  is  a  bifurcation  parameter.  For  small  perturbations  or 
noise,  the  dynamical  behavior  of  the  system  can  be  described  by  the  linearized 
equation  near  the  equilibrium  point  Xq.  The  linearized  system  corresponding  to 
(3.5)  with  a  small  noise  forcing  N(t)  is  given  by 

x  =  Df(x o,  /.i)x  +  N(t.)  (3.6) 

where  x  =  x  —  x0  and  N(t )  E  R"  is  a  vector  white  Gaussian  noise  having  zero 
mean.  For  the  results  of  the  linearized  analysis  to  have  any  bearing  on  the 
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power  spectrum  power  spectrum 


(a)saddle-node  or  transcritical 


frequency 
(c)period  doubling 


frequency 


(b)pitchfork 


Figure  3.1:  Power  spectrum  for  system  operating  at  a  limit  cycle  near  bifurcation 
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original  nonlinear  model,  we  must  assume  that  the  noise  is  of  small  amplitude. 
This  assumption  of  small  noise  will  be  made  clear  below,  in  terms  of  smallness 
of  correlation  and  cross-correlation  coefficients. 

The  noise  N(t)  can  occur  naturally  or  can  be  injected  using  available  controls. 
To  facilitate  consideration  of  cases  in  which  the  noise  is  intentionally  injected, 
we  write  N(t)  in  the  general  form 

N{t)  =  Bn(t)  (3.7) 

where  B  e  i?nxm  and  n(t)  €  Rm  is  a  vector  white  Gaussian  noise.  This  notation 
allows  us  to  easily  consider  cases  in  which  noise  is  injected  in  different  equations 
through  available  actuation  means.  The  noise  vector  N(t)  has  zero  mean  as  long 
as  n(t )  has  zero  mean. 

We  view  the  system  (3.6)  as  being  in  steady  state  and  driven  only  by  the 
noise  process.  Thus,  we  solve  for  the  evolution  of  the  state  assuming  a  zero 
initial  condition.  The  solution  of  equation  (3.6)  with  a  zero  initial  condition  is 

x{t)  =  eAt  f  e~AsN(s)cls  (3.8) 

Jo 

where  A  :=  DF(x o).  For  our  analysis,  we  assume  that  xo  is  an  asymptotically 
stable  equilibrium  point,  i.e.,  all  the  eigenvalues  of  A  have  negative  real  part. 
We  can  express  (3.8)  in  terms  of  the  eigenvectors  and  eigenvalues  of  A: 

n  rt  n 

x(t)  =  y~)  P  /  ^  lk N (s)eXkS rk dseXjtF 
j= i  Jo  k= 1 

/V  =  5tj 

where  /•'  and  /'  are  right  and  left  eigenvectors  of  A  corresponding  to  eigenvalue 
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A j,  respectively  and  where  5^  is  the  Ivronecker  delta, 

0  i  ±  j 
1  i  =  j 

Thus,  the  i-th  component  of  x(t)  is  given  by 


Sij  — 


(3.9) 


n  rt  n 

Xi(t)  =  EeXjtr!  Jo  E  lkN{s)e~XkSrkds 


(3.10) 

3= 1  "u 

Since  the  power  spectrum  is  the  Fourier  transformation  of  the  auto-covariance 
function,  we  calculate  auto-covariance  for  .r,,:  ( / ) . 


(xi{t)xi{t  +  t))  =  EE 


e\jte\k(t+T)  r3rk 


k  f*+T  f* 


e~\jSl  e~\kS2 


j=  1  fc  =  l 
n  n 

E  E  ^(^0(Sl)^(.s'2))d.Slds2 

o=l  p=l 


Let  the  noise  have  correlation  characteristic 


(Ni(t)Nj(t  +  t))  =  i \jS(T) 


(3,11) 


where  h(-)  is  the  Dirac  delta  function  and  the  ry,-  are  constants  for  all  i .  j .  More¬ 
over,  v,j.j  should  be  small  enough  such  that  linearized  analysis  is  valid.  Again, 
(3.11)  is  satisfied  for  the  linear  affine  input  system  as  long  as  n(t)  satisfies 


■  ••  •  +  •  =  •  ••  • 


(3,12) 


where  the  (/A  are  constants  for  all  i,j.  Then 

n  n 

(.r,- (/)./';(/  I  -))  =  Y.Y.(XA?XiA,':) 


k  f*+T  f* 


e_Aj'Sle_A*S2 


j= 1  k=l 
n  n 

EE  lj/pVoA8l  -  S2)clsids2 

0=1  p=  1 


ft 

=  EE  (  \jt(  Aj.S/  •  r)r.i f.k  /  e-\jSe~\ks 

j=lk=l  1 

n  n 

■  ZEfr** 

0=1  p-  I 


(3.13) 
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For  a  dynamic  system  which  depends  on  a  single  parameter,  there  are  two 
types  of  typical  bifurcation  from  a  nominal  equilibrium  point.  One  is  stationary 
bifurcation  in  which  case  a  new  equilibrium  emerges  or  the  original  equilibrium 
point  suddenly  disappears  at  the  bifurcation.  The  other  is  Hopf  bifurcation, 
where  a  periodic  orbit  emerges  from  the  equilibrium  point  at  bifurcation.  For 
stationary  bifurcation,  a  real  eigenvalue  of  the  linearized  system  becomes  zero 
as  the  parameter  varies.  For  Hopf  bifurcation,  a  complex  conjugate  pair  of 
eigenvalues  crosses  the  imaginary  axis. 

Consider  the  Hopf  bifurcation  case  first.  Assume  that  a  complex  conjugate 
pair  of  eigenvalues  (denote  them  as  A  =  Ai,  A  =  A2)  close  to  the  imaginary  axis 
has  relatively  smaller  negative  real  part  in  absolute  value  compared  to  other 
system  eigenvalues: 

|i2e(A1)|,|i2e(A2)|  <  |i?e(A.t)|  (3.14) 


for  i  =  3 . //.  Since  the  integrand  in  (3.13)  is  the  product  of  decaying  ex¬ 

ponentials  (due  to  the  asymptotic  stability  assumption)  and  a  bounded  value, 
terms  involving  Ai  and  A2  dominate  (3.13)  for  large  t : 


(Xi{t)Xi{t  +  t)) 


eAi(2  t+r) 


0\2  (2  t+r) 


If 

J  0  4  —  1  U—  1 

(n2)2  f 


j= 1 k= 1 

l  n  n 

— 2A2  s 


£  £ 

j- 1  k=l 


/l  n  n 

(Al  +A2)5 


0\2  (t+T )  +  A; 


'  r2  [ 

%  %  L 

Jo 

:'r)rv  f 

1  1  Jo 


E  E  fy&ijds 

j= 1  k= 1 


l  n  n 

e-a,+«.  E  E  iulvijds 


4=1  k=l 

The  power  spectrum  is  measured  with  the  use  of  a  spectrum  analyzer  and 
most  practical  spectrum  analyzers  perform  both  an  ensemble  average  and  a  time 
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average.  Thus,  the  final  auto-covariance  function  is 

Cu{t)  =  Re[(xi{t)xi{t  +  t))*]  (3.15) 

~  2e“eT  e-fT(ecos(a;r) -a>sin(a>r)) 

=  =■—  cosier) +  T-[ - - ]  (3.16) 

where  <  ■  >'  denotes  the  mean  over  time  t,  and  we  replaced  Ai  =  —  e  +  ico  and 
A2  =  —  e  —  tio  with  e,u>  >  0.  Also,  H  and  T  are 


:= 


i=l  k= 1 


T  :=  EE?fe(^)2  +  EEfe(^)2 

j=l  A:=l  j=l  A,:=l 

Finally,  taking  the  Fourier  transform  of  equation  (3.16)  yields  the  desired 
power  spectrum: 


n{V)  1-1  {jr/  +  e)2  +  uj2 

+  T[ _ e(^  +  f) _ I _ 13.17) 

(e2  +  uj2){{jr]  +  e)2  +  u2))  (e2  +  uJ2)((ji]  +  e)2  +  u;2))* 

The  magnitude  of  Su(r])  is  maximum  at  rj  =  lo  and  the  maximum  grows  without 

bound  as  e  — >■  0.  Moreover,  as  noise  power  (as  measured  by  the  )  increases,  the 

magnitude  of  Su(r])  also  increases.  However,  since  2  and  T  affect  Su(r])  linearly 

and  uniformly  over  frequency  rj,  the  shape  of  the  magnitude  S'**  (77)  doesn’t  change 

with  increasing  different  noise  power.  Of  course,  we  have  assumed  that  the  noise 

is  of  small  amplitude,  so  we  cannot  actually  allow  the  zy?-  to  increase  without 

bound. 

Fig.  3.2  shows  the  magnitude  of  Su(rj)  for  u  =  10.  Note  the  sharp  peak 
around  lo  =  10  that  appears  as  e  — >■  0.  From  this  observation,  we  can  conclude 
that  the  power  spectrum  peak  near  the  bifurcation  located  at  to ,  and  the  mag¬ 
nitude  of  this  peak  grows  as  e  approaches  to  zero.  This  property  will  be  used  as 
a  precursor  signaling  the  closeness  to  Hopf  bifurcation. 
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To  study  the  impact  of  noise  near  a  stationary  bifurcation,  assume  that 
a  real  eigenvalue  occurs  close  to  zero  (denote  it  as  A  =  Ai)  and  that  it  has 
relatively  smaller  negative  real  part  in  absolute  value  compared  to  the  other 
system  eigenvalues: 

| Ail  <C  |-Re(A,;)|  (3.18) 

for  i  =  2, . . . ,  n.  Due  to  (3.18),  terms  with  j  =  1  and  A:  =  1  dominate  the 
expression  (3.13)  for  large  t ,  so  that 

n  n 

£  £  l)l>nds 

j  I  k= 1 

Taking  the  time  average,  we  get  auto-covariance  function 

Cu{t)  :=  (xi(t)xi(t  +  r))t  (3.19) 

=  EDW(rJ)JT  (3-20) 

j= 1  k= 1 

Fourier  transformation  of  (3.20)  gives  the  desired  power  spectrum: 

n  n  i 

SU(V)  =  [E  E  l^r^2e{e  +  Jv)  (3'21) 

This  equation  shows  that  the  magnitude  of  the  power  spectrum  peak  grows  as 
e  approaches  to  zero  and  the  location  of  this  peak  is  r/  =  0.  Fig.  3.3  shows  the 
magnitude  of  Su(rj)  (3.21).  Note  the  sharp  growing  peak  around  uj  =  0  as  e  — )■  0. 


(Xi{t)Xi{t  +  t)) 


eX  i(2*+T)(rl)2 
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Chapter  4 


Monitoring  System  for  Detection  of  Proximity 
to  Stationary  and  Hopf  Bifurcation 


As  seen  in  the  proceeding  chapter,  we  can  expect  to  observe  a  growing  peak 
in  the  power  spectrum  of  a  measured  output  of  a  nonlinear  system  with  white 
Gaussian  noise  input  as  the  system  approaches  a  bifurcation.  In  the  case  of 
Hopf  bifurcation,  the  location  of  the  power  spectrum  peak  coincides  with  the 
imaginary  axis  crossing  frequency  of  the  critical  eigenvalues.  In  the  case  of 
stationary  bifurcation  the  power  spectrum  peak  occurs  at  zero  frequency.  In 
this  chapter,  we  start  with  these  observations  and  develop  monitoring  system 
that  yields  useful  precursors  if  proximity  to  bifurcation. 

4.1  Transforming  Stationary  Bifurcation  to  Hopf 
Bifurcation 

In  this  section,  introduce  a  the  monitoring  system  that  ,  when  augmented  to  a 
system  undergoing  stationary  bifurcation,  results  in  a  new  system  that  undergoes 
Hopf  bifurcation.  Hopf  bifurcation  is  easier  to  detect  than  stationary  bifurcation 
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with  the  and  of  noisy  precursors.  Hopf  bifurcation  shows  a  power  spectrum  peak 
at  a  non  zero  frequency. 

Consider  the  system 

&*=f(x,n)  (4.1) 

and  suppose  following  assumptions  hold. 

(Al)  The  origin  is  an  equilibrium  point  of  system  (4.1)  for  all  values  of  //,. 

(A2)  System  (4.1)  undergoes  stationary  bifurcation  at  //  =  jic.  (i.e.  there  is  a 
simple  eigenvalue  A(/i)  of  such  that  for  some  value  y  =  jic. 

K»c)  =  0  and  ^dT  +  °) 

(A3)  All  other  eigenvalues  of  Df( 0,nc)  are  in  the  open  left  half  complex  plane. 
Let  the  augmented  system  corresponding  to  (4.1)  be 


Xi  = 

f{x,  h)  ~  chi 

Vi  = 

CXi 

(4,2) 

Here,  y  E  Rn ,  c  E  R  and  i  =  1,2,... 

,  n. 

Proposition  4.1  Under  assumptions  (A1)-(A3),  the  augmented  system  (4- 2) 
experiences  Hopf  bifurcation  at  fi  =  fic.  Moreover,  if  for  any  value  of  fi  the 
origin  of  the  original  system  (4-1)  is  asymptotically  stable  (resp.  unstable),  then 
the  origin  is  asymptotically  stable  (resp.  unstable)  for  the  augmented  system 
(12). 

Proof:  Denote  by  A  the  Jacobian  matrix  for  system  (4.1)  at  the  origin.  Clearly, 
the  origin  (0,0)  in  R2n  is  an  equilibrium  point  of  the  augmented  system  (4.2). 
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The  Jacobian  matrix  of  the  augmented  system  (4.2)  at  the  origin  is 

A  -cl 

D  =  (4.3) 

cl  0 

Let  a  be  any  eigenvalue  of  A  and  r  the  corresponding  right  eigenvector.  Also, 
denote  by  A  any  eigenvalue  of  D  and  the  associated  eigenvector  by  v  =  [tq  'y2]T. 
Then, 


Aui 

=  Ac  1  —  CV  2 

(4.4) 

to 

=  CV  1 

(4.5) 

We  seek  a  solution  for  which  V\  = 

From  (4.5),  we  have 

c 

C2  =  ^ 

(4.6) 

Substituting  (4.6)  into  (4.4)  and  using  r  ^  0,  we  get 

A2 

—  a\  (?  —  0 

(4.7) 

Thus,  any  eigenvalue  a  of  A  has  corresponding  to  it  two  eigenvalues  of  D,  which 
are  the  solutions  of  quadratic  equation,  above: 


Thus,  the  eigenvalues  of  the  that  Jacobian  matrix  of  the  augmented  system  (4.2) 
are 

(\j  ±  J —  4c2 

'W  1.2;  - 2 -  *  =  l,2,...,n  (4.9) 

where  cq,  i  =  1, . . . ,  n  are  the  eigenvalues  of  A.  Let  the  eigenvalue  of  A  that 
becomes  0  be  ot\.  At  ^  =  /%  the  eigenvalues  of  the  augmented  system  associated 
with  a.i  are  (using  (4.7))  a  pair  of  pure  imaginary  eigenvalues  at  /ic 

Ai,A  2  =  ±cj  (4.10) 
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Note  that  the  pair  of  pure  imaginary  eigenvalues  (4.10)  depends  on  c. 

For  a  Hopf  bifurcation  to  occur,  the  crossing  complex  conjugate  pair  of  eigen¬ 
values  should  satisfy  the  transversality  condition.  From  the  quadratic  equation 
(4.7)  and  using  the  fact  that  ot\  =  0  at  //  =  (ic,  we  have 

dRe(Xi)  1  da\ 

dji  2  dji  '  y 

At  ^  =  iic.  Since  a.\  =  0  and  ^  ^  0  at  //  =  jir  from  assumption  (A2),  (4.11) 
implies  dRe}^  =  7^  0  (i.e.,  the  transversality  condition  hold  for  system 

(4.2)).  Therefore,  the  augmented  system  (4.2)  undergoes  Hopf  bifurcation  from 
the  origin  at  //  =  (ic. 

The  last  step  in  the  proof  consists  in  showing  that  all  other  eigenvalues  of  the 
matrix  D  are  in  the  open  left  half  complex  plane.  Any  pair  of  eigenvalues  of  D 
can  be  obtained  from  (4.9).  For  a  real  eigenvalue  of  A ,  it  is  clear  from  (4.9)  that 
the  corresponding  pair  of  eigenvalues  of  D  have  negative  real  part  if  cq  <  0  since 
cq  <  Re[  Jaf  —  4c2].  For  a  complex  conjugate  pair  of  eigenvalues  of  A  (denote 


them  as  7,7),  we  have  following  the  two  equations: 

A2  — 7A  +  C2  =  0  (4.12) 

A2  — 7A  +  C2  =  0  (4.13) 

Multiply  (4.12)  and  (4.13)  to  get  the  following  fourth  order  equation: 

A4  —  (7  +  7)  A3  +  (2c2  +  77)A2  —  c2(7  +  j)X  +  c4  =  0  (4-14) 

Denoting  7  =  a  +  bj  and  7  =  a  —  bj ,  equation  (4.14)  simplifies  to 

A4  -  2a A3  +  (2c2  +  a2  +  b2) A2  -  2ac2A  +  c4  =  0  (4.15) 
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Applying  the  Routh-Hurwitz  criterion  [25]  to  (4.15),  we  obtain  the  Routh  array 


s4 

1 

2c2  +  a2  +  b2 

e4 

.s3 

—  2a 

—2ac2 

0 

s2 

3c2  +  a2  +  b 2 

c4 

0 

s1 

—  2ac2(2c2+a62+b2) 
3c2+a2+62 

0 

0 

s° 

c4 

0 

0 

From  assumption  (A3),  a  <  0.  This  guarantees  that  all  the  entries  in  the  first 
column  of  Routh  array  are  positive.  Therefore,  all  eigenvalues  of  the  Jacobian 
matrix  of  the  augmented  system  have  negative  real  part. 

From  forging  discussion,  it  is  also  clear  that  any  eigenvalue  of  A  has  positive 
real  part,  then  the  corresponding  eigenvalues  of  D  also  have  positive  real  part. 
This  proves  that  if  the  original  system  is  unstable,  then  the  augmented  system 
is  also  unstable.  □ 

Note  that  by  changing  the  value  c  in  equation  (4.2),  we  can  control  the 
crossing  frequency  of  the  complex  conjugate  pair  of  eigenvalues  of  the  augmented 
system.  Thus,  for  detecting  stationary  bifurcation,  we  only  need  to  monitor  a 
frequency  band  around  the  value  c  ,  which  we  control. 

There  are  some  other  advantages  of  our  monitoring  system.  The  augmented 
system  (4.2)  has  the  same  critical  parameter  value  (fic)  as  the  original  system. 
This  is  actually  not  a  luxury  but  a  necessity  for  the  system  to  be  practically 
useful.  In  addition,  the  final  part  of  the  proof  shows  that  augmenting  the  state 
y,  and  applying  the  feedback  q/j  to  the  original  system  does  not  change  stability 
of  the  system.  Moreover,  use  of  design  does  not  require  knowledge  of  the  original 
system.  However,  there  are  some  restriction  to  apply  to  suggested  system  to 
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general  nonlinear  system.  We  will  discuss  this  in  the  some  of  the  following- 
sections. 

We  assumed  that  in  the  original  system  the  stationary  bifurcation,  the  transver- 
sality  condition  is  satisfied.  It  is  not  the  case  for  saddle  node  bifurcation.  For  a 
saddle-node  bifurcation,  the  augmented  system  of  this  section  results  in  a  degen¬ 
erate  Hopf  bifurcation.  The  possible  bifurcation  diagrams  for  degenerate  Hopf 
bifurcation  are  more  complex  than  for  Hopf  bifurcation  [47],  [32],  However,  for 
the  purposes  of  the  power  spectrum  precursor  it  is  enough  that  a  new  peak  in 
the  power  spectrum  appears  near  bifurcation. 

4.2  Detection  Hopf  Bifurcation  using  Monitor¬ 
ing  System 

In  this  section,  we  consider  the  effect  of  the  monitoring  system  of  the  proceed¬ 
ing  section  on  a  system  which  undergoes  Hopf  bifurcation  instead  of  stationary 
bifurcation.  Consider  again  the  system  (4.1),  repeated  here  for  convenience: 

■r-  fix.  ii)  (4.16) 

(Al)7  The  origin  is  an  equilibrium  point  of  (4.16)  for  all  values  of  //,. 

(A2)  System  (4.16)  undergoes  a  Hopf  bifurcation  from  the  origin  for  //  =  /ic, 

(A3)'  All  other  eigenvalues  of  Df{ 0,^c)  are  in  the  open  left  half  complex  plane. 
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Let  the  augmented  system  be 


Xi  =  fi(X./l)  <  'I/ ; 

Vi  =  cXi  (4.17) 

where  x  G  Rn ,  y  G  Rn ,  c  G  R,  and  i  =  1,2,...,  n. 

Proposition  4.2  Uiider  the  assumptions  (Al)' -(AS)' ,  the  augmented  system 
(4-17)  undergoes  a  codimension  two  bifurcation  at  y  =  yc,  in  which  two  complex 
conjugate  pair  of  eigenvalues  cross  the  imaginary  axis.  Moreover,  for  any  value 
of  y  if  the  origin  of  the  original  system  is  asymptotically  stable  (resp.  unsta¬ 
ble),  then  the  origin  is  asymptotically  stable  (resp.  unstable)  for  the  augmented 
system. 

Proof:  First,  we  show  that  the  augmented  system  has  two  pair  of  pure  imaginary 
eigenvalues  at  the  origin  for  //,  =  //,c,  and  that  these  eigenvalues  satisfy  the 
transversality  condition. 

From  assumption  (A2)\  the  Jacobian  matrix  of  the  original  system  at  the  ori¬ 
gin  has  a  pair  of  pure  imaginary  conjugate  eigenvalues  (denote  them  by  jco.  —juf) 
for  pL  =  fie.  From  the  proof  of  Proposition  4.1,  it  is  clear  that  each  of  these  eigen¬ 
values  results  in  a  pair  of  eigenvalues  for  the  augmented  system  which  are  the 
solutions  of  the  following  equations 

X2-jcv\  +  e2  =  0  (4.18) 

A2  +  ju)X  +  c2  =  0  (4.19) 

By  multiplying  the  equations  above,  we  get  a  fourth  order  equation  the  solutions 
of  which  are  eigenvalues  of  augmented  system: 

A4  +  (2c2  +  uj2)X2  +  c4  =  0  (4.20) 
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The  four  solutions  of  the  equation  above  are  given  by 


A  _  ±  \/~2c2  -u2±  \/rTc].c’  +  CJ4 

V2 

Note  that  2c2  +  co2  >  \J 4 c2uj2  +  uja  for  all  o,  u  G  R.  Therefore,  the  Jacobian 
matrix  of  the  augmented  system  has  two  pairs  of  pure  imaginary  conjugate 
eigenvalues  at  the  critical  parameter  value. 

To  check  the  transversality  condition,  consider  the  eigenvalues  for  //  near  gc. 
Near  g  =  g,c,  we  have  the  following  fourth  order  equation  the  solutions  of  which 
result  from  the  pair  of  complex  conjugate  eigenvalues  (a  ±  toj)  of  the  original 
system  (see  (4.14)): 

A4  -  2aA3  +  (2c2  +  a2  +  uj2) A2  -  2c2  a  +  e4  =  0  (4.22) 


Since  //  is  close  to  /ic,  by  continuity  it  follows  from  above  that  this  equation  has 
two  pairs  of  complex  conjugate  eigenvalues  as  its  solutions  for  ^  near  /ic.  Denote 
these  as  e±fj,g±hj.  Then,  we  have  (4.24)  from  relationship  between  solutions 
of  equations  and  its  coefficients. 

4 

-2a  =  A i  =  e  +  g  (4.23) 

2=1 
4 

—  2c2 a  =  ^  \iXjXk  =  2 g(e2  +  f2)  +  2 e(g2  +  h 2)  (4.24) 

i,j,k=l 

i+j+k 

where  A,.  X}.  A*:  are  roots  of  (4.22).  Take  the  derivative  of  both  sides  of  the 
equations  above  with  respect  to  ji  and  evaluate  at  g  =  jic  (e  =  g  =  0,  f3  = 
0,7 2  =  c  -  a),  giving 


de 


dg 


dq 

+  7 

j2  dg 
dg 


^ da 
dg 
2  da 
dg 


(4.26) 
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Solve  equations  (4.25)  for  Also,  h  and  /  are  not  0  at  the  critical  point 

from  (4.21)  and  /  ^  h  at  the  critical  point  if  c  ^  0.  These  conditions  guarantee 
that  if  g  ^  0  ,  then  g  ^  0,  J  ^  0. 

The  last  step  in  the  proof  consists  in  showing  that  all  other  eigenvalues  of 
Jacobian  matrix  D  of  (4.17)  in  the  open  left  half  complex  plane.  Note  that  we 
have  same  form  of  matrix  D  as  in  proof  of  Proposition  4.1. 


D 


A  -cl 
cl  0 


(4.26) 


where  all  noncritical  eigenvalues  of  A  have  negative  real  part.  We  can  use  the 
same  procedure  as  in  Proposition  4.1  to  prove  that  if  all  noncritical  eigenvalues  of 
A  have  negative  real  part,  then  all  corresponding  eigenvalues  of  D  have  negative 
real  part. 

It  is  also  clear  that  if  any  eigenvalue  of  A  has  positive  real  part,  then  the 
corresponding  eigenvalues  of  D  also  have  positive  real  part.  This  proves  that  if 
the  original  system  is  unstable,  then  the  augmented  system  is  also  unstable.  □ 


Since  two  pairs  of  eigenvalues  of  the  augmented  system  cross  the  imaginary 
axis  at  the  critical  parameter  value,  we  can  expect  to  see  two  peaks  in  the  power 
spectrum  as  the  system  nears  the  bifurcation  point.  From  (4.21),  we  see  that 
values  of  the  pairs  of  imaginary  eigenvalues  at  criticality  depend  on  c.  Hence,  we 
can  change  the  location  of  the  power  spectrum  peaks  by  changing  c.  Moreover, 
we  can  predict  the  exact  locations  of  the  peaks  if  the  imaginary  axis  crossing- 
pair  of  eigenvalues  of  the  original  system  are  known. 

However,  the  augmented  system  undergoes  a  codimension  two  bifurcation. 
Thus,  the  post-bifurcation  behavior  is  complicated.  The  nature  of  the  bifurcation 
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depends  strongly  on  The  detailed  possible  bifurcation  diagrams  for  the 

corresponding  normal  form  can  be  found  in  [13],  [27],  [38]  and  [31]. 

4.3  Reduced  Order  Monitoring  System 

There  are  some  drawbacks  of  the  monitoring  system  approach  presented  in  the 
preceding  two  sections.  Most  crucial  of  these  is  that  it  requires  that  all  of  the 
original  system  states  be  measured  and  fed  back  through  the  augmented  system. 
This  is  not  be  practical  for  many  systems.  In  this  section,  we  show  that  for 
some  classes  of  nonlinear  systems  we  can  apply  a  similar  monitoring  system 
without  needing  feedback  of  all  the  original  system  states.  This  is  achieved  by 
supposing  the  original  system  has  an  inherent  time-scale  structure,  i.e.,  that 
it  is  in  singularly  perturbed  form  or  the  system  with  affine  input  has  linearly 
decoupled  Jacobian  matrix. 

4.3.1  Monitoring  system  for  singularly  perturbed  system 

Consider  a  general  singular  perturbation  problem  first.  Let  the  full  system  be 
of  the  form 


x  =  f(x.z,/i.<) 

ez  =  g(x,  z,  g,  e)  (4-27) 

where  x  G  Rn ,  z  G  Rm ,  //.  e  G  R  and  e  is  small  but  positive.  The  reduced  system 
is  obtained  by  formally  setting  e  =  0  in  (4.27),  giving 

x  =  f(x,  z,  fi,  0) 

0  u(x.z.iiA))  (4.28) 
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Let  ?77-o  =  (0,  -2q)  be  an  equilibrium  point  of  the  reduced  system.  Also,  assume 


(HI)  777o  =  (0,^o)  is  an  equilibrium  point  of  (4.28)  for  all  values  of  //. 

(H2)  /,</,  are  Cr(r  >  5)  in  x,  z,  //,,  e  in  a  neighborhood  of  (m.0,  0,  0). 

(H3)  det(JD2^(0,-20,0,0))  ^  0. 

(H4)  The  reduced  system  undergoes  a  stationary  bifurcation  at  mo  for  the  crit¬ 
ical  parameter  value  fi  =  \ic. 

(H5)  No  eigenvalue  of  Dzg(0,  z0, 0,  0)  has  zero  real  part. 

Let  the  augmented  system  corresponding  to  (4.27)  be 

x  =  f(x,  zrn,  e)  —  cy 
y  =  ex 

ez  =  g(x,z,n,e)  (4.29) 

Proposition  4.3  Let  (H1)-(H5)  above  hold.  Then  there  is  an  e o  >  0  and  for 
each  e  G  [0,e0]  the  augmented  system  (f.29)  undergoes  a  H op f  bifurcation  at  an 
equilibrium  m-o°’e  'near  m o  for  a  critical  parameter  value  ffc  near  gc. 

Proof:  By  virtue  of  Theorem  2.4  an  persistence  of  Hopf  bifurcation  under  singular 
perturbation,  we  only  need  to  show  that  the  reduced  system  corresponding  to 
(4.29)  undergoes  a  Hopf  bifurcation  at  (0,0,  z0)  at  the  critical  parameter  value 
//  =  gc.  The  reduced  system 

x  =  f(x,z,/i,0)  -  cy 
y  =  cx 

0  =  g(x,z,n,  0)  (4.30) 
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Since  the  original  reduced  system  (4.28)  undergoes  a  stationary  bifurcation,  we 
can  apply  Proposition  4.1  to  (4.30)  to  show  that  the  reduced  augmented  system 
(4.30)  undergoes  Hopf  bifurcation  at  the  critical  parameter  value  //  =  fic.  And 
then,  The  remainder  of  proof  follows  as  for  the  Theorem  2.4. □ 

Proposition  4.3  is  useful  because  it  implies  that  we  only  have  to  augment 
and  feed  back  slow  states  in  a  two-time  scale  system  to  transform  stationary 
bifurcation  into  Hopf  bifurcation. 

4.3.2  Monitoring  system  for  the  system  with  linearly  de¬ 
coupled  Jacobian  matrix 

In  this  section,  we  consider  a  nonlinear  affine  control  system  which  has  linearly 
decoupled  Jacobian  matrix  at  origin.  For  those  system,  we  introduce  new  re¬ 
duced  order  monitoring  system. 

Let  this  nonlinear  affine  control  system  be 

n 

x  =  f(x,  n)  +  9iix)ui  (4.31) 

1=1 

where  is  smooth  function  of  Xi  and  each  rq  is  an  input  which  will  be  used 
to  feedback  the  augmented  states.  If  the  controls  iq  =  0  for  i  =  1, . . , ..in.  the 
system  becomes 

x  =  f{x,n)  (4.32) 

The  following  assumptions  on  (4.31)  and  (4.32)  are  used  in  Proposition  4.4  below. 
(Bl)  The  origin  is  an  equilibrium  point  of  (4.32)  for  all  //,. 

(B2)  System  (4.32)  undergoes  stationary  bifurcation  at  n  =  jic  (i.e.,  the  uncon¬ 
trolled  system  experiences  stationary  bifurcation). 
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(B3)  Besides  the  critical  zero  eigenvalue,  all  other  eigenvalues  of  the  system 
linearization  are  in  the  open  left  half  complex  plane  for  p  =  p.c. 


(B4)  The  Jacobian  matrix  of  the  system  at  origin  for  //  =  pc  has  the  form 


A  :  = 


An  0 
A21  A22 


(4.33) 


where  An  G  Rlxl  with  l  <  m.  Moreover,  the  zero  eigenvalue  of  the  system 
at  criticality  is  an  eigenvalue  of  the  An. 


(B5)  The  subspace  D(x )  =  span{gfix)  :  i  G  1 has  constant  dimension 
k  >  l  for  all  x  in  a  neighborhood  of  the  origin. 


Proposition  4.4  Introduce  the  augmented  system 

n 

x  =  f{x,g)  +  Y,9i{x)u*. 


y  =  c 


I  0 


2=1 


X 


(4.34) 


where  x  G  Rn ,  y  G  Rl ,  c  G  R.  and  I  is  l  x  l  identity  matrix.  Under  assumptions 
(B1)-(B5),  the  augmented  system  (4-34)  undergoes  Hopf  bifurcation  at  p  =  pc 
with  input  u*  =  —cKy.  Here,  c  G  R  and  I\  G  Rmxl  satisfies 


G{0)I< 


I 

0 


(4.35) 


Here,  G(0)  =  [^(O),  •  •  • ,  <?OT(0)]  and  I  is  l  x  l  identity  matrix.  In  addition,  if 
origin  is  asymptotically  stable  for  the  original  system,  then  the  origin  is  asymp¬ 
totically  stable  for  the  augmented  system. 
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Proof:  Because  of  assumption  (B4),  there  exists  an  K  satisfies  (4.35).  Then,  the 
Jacobian  matrix  of  the  augmented  system  with  u*  is 

An  0  —cl 

D=  A21  H22  0  (4.36) 

cl  0  0 

Since  interchanging  rows  and  columns  does  not  effect  the  eigenvalues  of  a  matrix, 
we  reformulate  D  as  follows: 

Di  0 

D  =  (4.37) 

A21  0  A22 

where 

An  -cl 
cl  0 

Since  the  set  of  eigenvalues  of  D  is  the  union  of  the  eigenvalues  of  Di  and  of 
7l22,  we  can  now  show  that  the  system  undergoes  Hopf  bifurcation  by  applying 
Proposition  4.1  to  Di.  It  is  also  clear  that  the  origin  is  asymptotically  stable 
for  the  augmented  system  as  long  as  it  is  asymptotically  stable  for  the  original 
system  from  the  proof  of  Proposition  4.1.  □ 

If  we  have  system  which  is  the  case  with  (4.31)  under  assumption  (B4)  and  if 
system  experiences  stationary  bifurcation,  then  we  can  transform  the  bifurcation 
into  a  Hopf  bifurcation  with  reduced  order  augmented  states.  In  the  extreme 
case,  we  need  only  a  single  augmented  state  if  the  dimension  of  An  is  one.  It  also 
have  advantage  over  multi-dimension  An  which  will  be  illustrated  with  some  ex¬ 
ample  in  section  5.2.  For  those  system  which  is  not  decoupled  as  (B4),  we  could 
try  coordinate  transform  which  render  the  system  into  desired  form.  However, 
detailed  knowledge  of  the  system  is  needed  to  find  such  a  transformation. 


(4.38) 


53 


4.4  Stability  of  Bifurcated  Periodic  Solution  for 


the  Augmented  System 

In  this  section,  we  consider  the  stability  of  bifurcated  solutions  of  the  augmented 
system.  We  have  shown  in  previous  section  that  the  proposed  monitoring  system 
renders  stationary  bifurcation  into  Hopf  bifurcation. 

We  will  investigate  in  this  section  the  stability  of  the  bifurcated  periodic- 
orbit  if  the  original  system  undergoes  either  supercritical  or  subcritical  pitchfork 
bifurcation. 

Let  the  system  be 

x  =  f(x,fi)  (4.39) 

where  x  G  Rn  is  the  state  vector  and  fi  G  R  is  the  bifurcation  parameter. 
Suppose  that  at  the  critical  parameter  value  /./.  =  fi-c,  the  Jacobian  matrix  of 
(4.39)  evaluated  at  the  equilibrium  point  xo  =  0  has  one  zero  eigenvalue. 

Consider  the  case  in  which  n  =  1,  that  is,  let  the  dimension  of  the  state  vector 
of  (4.39)  be  one.  Also,  suppose  system  (4.39)  undergoes  a  pitchfork  bifurcation. 
It  is  easy  to  see  that  left  (/)  and  right  (r)  eigenvector  corresponding  to  the  simple 
zero  eigenvalue  at  criticality  can  be  taken  as  any  nonzero  constants.  Set  r  =  1 
and  l  =  1  such  that  r  and  l  satisfy  Ir  =  1  and  the  first  component  of  r  is  one. 
Lemma  2.1  then  applies.  From  the  assumption  that  system  (4.39)  undergoes 
pitchfork  bifurcation,  we  have 

d2  f 

Pi  =  lQ(r ,  r )  =  -q^{ °)  =  °  (4.40) 

Thus,  Q(r,r)  is  zero.  Recall  that  /L  is  given  by 

P2  =  2l{2Q0(r1x2)  +  Coir,  r,  r)}  (4.41) 
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where  x,  y,  //,  G  R.  We  have  shown  that  the  augmented  system  undergoes  Hopf 
bifurcation  if  the  original  system  undergoes  pitchfork  bifurcation.  To  check  the 
stability  of  the  bifurcated  periodic  solution  of  (4.45),  we  only  have  to  check  the 
sign  of  the  Hopf  bifurcation  stability  coefficient  fa  (2.50).  At  criticality,  the 
Jacobian  matrix  of  (4.45)  is 

0  -c 

L0  :=  (4.46) 

c  0 

The  matrix  L0  has  an  eigenvalue  cj  with  corresponding  right  eigenvector  r  = 

r  iT  r  i 

1  _j  and  left  eigenvector  l  =  \  1  j  ,  respectively.  Eigenvalue  —  cj 

has  right  eigenvector  f  and  left  eigenvector  respectively.  Note  that  higher 
order  terms  only  come  from  f(x,  //)  and  they  are  not  a  function  of  y.  From  this 
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observation,  we  have 


Q((x,y),  (x,y)) 


where  we  used  the  equation  (4.40).  Equation  (4.47)  implies  that  Q(r,r)  and 
Q(r,r)  are  equal  to  zero.  Therefore,  a  =  0  and  6  =  0  are  solution  of  (2.48)  and 
(2.49),  respectively.  Hence,  13 2  of  (4.45)  becomes 

p2  =  2Re\)-lC(r,r,f)]  =  2--^(0)  (4.48) 

since  higher  order  terms  only  come  from  f(x,[i)  and  they  are  not  function  of  y. 

Note  that  a  sign  of  (4.44)  are  equal  to  a  sign  of  (4.48).  We  have  following- 
proposition  as  a  result. 

Proposition  4.5  Consider  the  system  (4-39)  is  of  first  order,  i.e.,  n=l.  If  the 
original  system  (4-39)  undergoes  a  supercritical  pitchfork  bifurcation  (respectively 
a  subcritical  pitchfork  bifurcation) ,  then  the  transformed  system  (4-45)  undergoes 
a  supercritical  Hopf  bifurcation  (respectively  a  subcritical  Hopf  bifurcation) . 

Next,  we  are  going  to  consider  the  case  when  n  >  2,  that  is,  dimension  of 
original  system  is  2  or  higher.  We  use  an  example  to  show  that  a  supercritical 
pitchfork  bifurcation  (resp.  a  subcritical  pitchfork  bifurcation)  need  not  result 
in  a  to  supercritical  Hopf  bifurcation  (resp.  a  subcritical  Hopf  bifurcation)  by 
using  the  monitoring  system  proposed  in  the  preceding  sections. 


1 

1 

Q5IQ5 

tor~h 

3 

0 

X 

2! 

x  y 

0 

0 

y 

(  1  <tL 


2  e,"  f° 

0 


(A 


v0/ 


(4.47) 
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To  this  end,  consider  the  example 


X\  =  —jlX  i  —  x\  +  X\X2 

X2  =  —X2  +  kx\  (4.49) 

where  //,  E  R  is  a  bifurcation  parameter  and  k  E  R  is  a  constant.  It  is  easy  to 
see  that  the  origin  is  an  equilibrium  point  for  all  parameter  values  //  and  that  a 
pitchfork  bifurcation  occurs  for  //  =  0.  Moreover,  a  simple  calculation  shows  that 
(3 2  for  this  pitchfork  bifurcation  is  -1.  This  implies  that  system  (4.49)  undergoes 
a  supercritical  pitchfork  bifurcation  at  //  =  0. 

The  augmented  system  corresponding  to  (4.49)  is 

x,\  =  —[ix  i  —  x\  +  X\X2  —  cyi 
{ji  =  ex  i 

x2  =  -X2  +  kx\  -  cy2 

2/2  =  CX2  (4.50) 

The  Jacobian  matrix  of  (4.50)  evaluated  at  the  origin  at  criticality  is 

0  -c  0  0 

c  0  0  0 

(4.51) 

0  0  -1  -c 

0  0  c  0 

This  matrix  has  eigenvalue  ci  with  corresponding  right  eigenvector  r  =  1  —  j  0  0 

and  left  eigenvector  l  =  \  1  j  0  0  •  The  eigenvalue  —  cj  has  right  eigen¬ 

vector  r  and  left  eigenvector  l.  The  Taylor  series  expansion  of  the  right  side  of 
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(4.50)  has  the  following  quadratic  and  cubic  terms: 


Q((x,y),(x,y))  = 


X\X2 

0 

hr2 

r\jJb  ^ 

0 


C({x,ti).(’x:y).(x.y))  = 


Therefore,  we  have 


Q(r,f)  =  Q(r,  r)  = 


—x 

0 

0 

0 

0 

0 

A: 

0 


C(r,  r,  r )  = 


-1 

0 

0 

0 


By  solving  (2.48)  and  (2.49),  we  have 


a  = 


b  = 


i  T 


0  0  0  -^ 


1  T 


0  0 


kj 


2j-3c  2j-3c 

Putting  these  values  into  (2.50),  (32  for  the  system  (4.50)  is 

A  =  -3  A; 


4  4  +  9c2 


(4.52) 


(4.53) 


(4.54) 


(4.55) 
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For  some  large  positive  k,  62  becomes  positive.  Thus,  the  augmented  system 
(4.50)  undergoes  a  subcritical  Hopf  bifurcation  for  such  value  of  k. 

Thus,  for  n  >  2  if  the  original  system  undergoes  a  supercritical  pitchfork 
bifurcation  (resp.  a  subcritical  pitchfork  bifurcation),  then  the  augmented  sys¬ 
tem  need  not  undergo  a  supercritical  Hopf  bifurcation  (resp.  a  subcritical  Hopf 
bifurcation). 

4.5  Nonlinear  Monitoring  System  Ensuring  Sta¬ 
bility 

In  this  section,  we  modify  our  monitoring  system  such  that  if  original  system 
undergoes  pitchfork  bifurcation  for  any  type  with  mild  assumption  which  will 
be  specified  later,  then  new  nonlinear  monitoring  system  undergoes  supercritical 
Hopf  bifurcation,  that  is,  a  bifurcated  periodic  solution  is  a  stable  limit  cycle. 

Consider  the  following  system  assumed  to  undergo  a  pitchfork  bifurcation  for 
k  kc- 

x  =  f(x,n)  (4.56) 

Here,  x  G  Rn  and  ^  G  R  is  the  bifurcation  parameter.  Denote  rs  and  ls  as  a 
right  eigenvector  and  left  eigenvector  corresponding  to  a  simple  zero  eigenvalue 
at  //  =  He,  respectively.  Take  first  component  of  rs  to  be  1  and  impose  the 
normalization  /srs  =  1. 

Consider  the  following  augmented  system 

Xi  =  fi{x,  k)  —  CXji 

Vi  =  cXi  -  mx\yi  (4.57) 
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where  m  is  a  real  constant.  At  criticality,  it  has  Jacobian  matrix 

A  -cl 

D  :=  (4.58) 

cl  0 

where  A  is  the  Jacobian  matrix  of  (4.56).  Employing  Proposition  4.1,  it  is  easy 
to  show  that  the  augmented  system  (4.57)  undergoes  Hopf  bifurcation  and  D 
has  eigenvalue  ±d.  Moreover,  the  right  and  left  eigenvectors  of  D  corresponding- 
eigenvalue  cj  are  given  by 

r  1 T 

r  =  rs  —jrs 

1  =  [l.  jl.]  («9) 

Also,  the  right  and  left  eigenvectors  corresponding  to  the  eigenvalue  —  cj  are 
given  by  r  and  l ,  respectively. 

The  stability  of  the  bifurcated  periodic  solution  of  augmented  system  is  de¬ 
termined  by  the  sign  of  fi-i  (2.50): 

P2  =  2Re[2lQ0(r\  a)  +  lQ0{f,  b)  +  ^/f(r.  r,  f)]  (4.60) 

Since  there  are  no  quadratic  terms  in  y,  in  the  augmented  system  (4.57),  /T  of 
(4.57)  simplifies  to 

P 2  =  2Re[2lsQf(r,a)  +  lsQf{f,b)  +  jlC(r,r,r)\  (4.61) 

where  Q/(-,  •)  denotes  a  quadratic  term  of  /(:r,yc).  Moreover,  we  can  simplify 

/C(r,  r,  f)  =  lsCf{r ,  r,  f)  +  jlsCy{r ,  r,  f ) 

771  . 

=  lsCf(r,r,  f)  -  —  J2iys{r])2  (4.62) 

z  i= l 

where  Cj(x ,  #,  .x)  are  cubic  part  of  f(x.  yc)  and  r*  and  l\  denotes  the  i-tli  compo¬ 
nent  of  rs  and  ls ,  respectively.  Since  the  first  component  of  r\  is  1  and  lsrs  =  1, 
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(4.62)  reduces  to 

rn 

lC(r ,  r,  f)  =  lsCf(r ,  r,  f)  — —  (4.63) 

Hence,  fa  becomes 

771 

A  =  2ife[2W,(r,a)  +  W/(f,6)  +  4C>(r,  r,f)]  -  -  (4.64) 

By  choosing  m  large  enough  and  positive,  we  can  make  fa  negative.  This  will 
imply  that  the  augmented  system  (4.57)  undergoes  a  supercritical  Hopf  bifurca¬ 
tion. 

Here,  we  suggested  only  one  of  many  possible  methods  which  renders  bifur¬ 
cated  periodic  orbit  stable.  Note  that  we  added  a  nonlinear  term  only  to  the 
dynamics  of  the  augmented  states  yt  not  the  physical  system  states  Xj.  This 
implies  that  we  have  freedom  to  choose  m  to  be  any  value. 

4.6  Monitoring  System  for  Equilibrium  Point 
Not  at  the  Origin 

We  have  shown  that  using  the  augmented  system  results  in  a  Hopf  bifurcation 
if  the  original  system  undergoes  a  stationary  bifurcation.  However,  the  scheme 
used  places  a  strict  requirement  on  the  system.  In  the  assumptions  of  Proposition 
4.1,  note  that  (Al)  demands  that  the  origin  be  an  equilibrium  point  of  the  system 
for  all  parameter  values.  This  assumption  is  invoked  so  that  the  equilibrium 
point  of  the  original  system  is  not  changed  through  adding  state  feedback.  In 
this  section,  we  attempt  to  remove  assumption  (Al).  It  will  become  clear  at 
the  end  of  the  section  that  removing  (Al)  comes  at  some  expense  in  terms  of 
simplicity  of  the  results. 
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Let’s  assume  that  equilibrium  point  of  system  (4.1)  of  interest  does  not  lie 
at  the  origin,  and  that  it  changes  as  //  varies.  Then  the  augmented  system  has 
equilibrium  point  that  differs  from  the  original  system.  As  a  matter  of  fact, 
(0,  ?/o)  becomes  a  new  equilibrium  point  for  augmented  system  (4.2),  where  1/0 
solves  /( 0,  //)  —  cy  =  0.  This  sudden  change  in  the  equilibrium  point  is  highly 
undesirable.  It  can  result  in  unstable  equilibrium  point,  or  in  some  cases  we  sim¬ 
ply  prefer  not  to  change  the  nominal  operating  point.  Therefore,  an  alternative 
augmented  system  for  monitoring  for  nearness  to  instability  is  presented  next. 
The  new  system  is  motivated  by  washout  filters,  discussed  in  Section  2.4. 

Consider  the  augmented  system 


Xi 

% 

0 

I 

'a 

II 

Vi 

=  CXi  +  a-Zi 

h 

=  Vi 

(4.65) 

where  i  =  1,  2, . . . ,  n  and  a,  c  £  R. 

Proposition  4.6  Assume  the  original  system  (4-1)  satisfies  (A2)  and  (A3)  with 
an  equilibrium  point  xo  not  necessarily  at  the  origin.  Then  the  augmented  system 
(4-65)  undergoes  a  codimension  2  bifurcation  at  p  =  pc.  At  criticality,  the 
linearization  of  (4-65)  possesses  one  simple  zero  eigenvalue  and  a  pair  of  pure 
imaginary  eigenvalues. 

Proof:  The  equilibrium  point  of  the  augmented  system  (4.65)  is  ( Xq ,  0,  Zo),  where 
zo  is  solution  of  cXi  +  az ,  =  0.  Note  that  new  augmented  system  keeps  xq  as  a 
component  of  this  equilibrium  point.  The  Jacobian  matrix  of  (4.65)  evaluated 
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at  this  equilibrium  point  is 


D 


A  —cl  0 
cl  0  al 


0 


I 


0 


(4.66) 


where  A  is  the  Jacobian  matrix  of  the  original  system  evaluated  at  xq.  Let  a 
be  any  eigenvalue  of  A  and  r  corresponding  eigenvector.  Also,  assume  A  is  an 
eigenvalue  of  D  with  eigenvector  v  =  [v{  v 2  vJ]T .  Then 


A  tii 

=  Av  1  —  CV  2 

(4.67) 

Xv2 

=  cv  i  +  av  3 

(4.68) 

Xv3 

=  v2 

(4.69) 

Attempt  a  solution  v  for  which  v\  =  r.  Solve  (4.68)  and  (4.69)  for  v2  and  v3  in 
terms  of  r,  we  get 


cA 


V2  = 


«3  = 


A2  —  a 


A2  —  a 

Substituting  the  equation  for  v2  into  (4.67),  gives 


A3  —  ccA  T  (c  —  a) A  +  cio,  —  0 


(4.70) 


Since  one  eigenvalue  of  A  becomes  0  at  //  =  /uc,  we  can  set  a  =  0  to  get  following- 
equation  for  the  expected  pair  of  eigenvalues  system  at  criticality: 

A3  +  (c2  -  a) A  =  0  (4.71) 

If  we  choose  a  <  0,  then  D  has  eigenvalues  0,  ±\/c2  —  a  j  which  correspond  to 
the  zero  eigenvalue  of  the  original  system  criticality. 
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Next,  we  check  the  transversality  condition.  Equation  (4.70)  which  corre¬ 
sponds  to  crossing  simple  real  eigenvalue  of  original  system  has  one  real  and  a 
pair  of  complex  conjugate  as  its  solution  near  the  critical  point.  Denote  S  as  the 
real  eigenvalue  and  /?  ±  7 j  as  the  pair  of  complex  conjugate  eigenvalue.  Solving 
this  notation  in  (4.70)  and  separating  real  and  imaginary  parts,  we  obtain 


5  +  2  *  f3  =  —a 

5{02  +  y2)  =  act  (4.72) 


Differentiating  these  equations  with  respect  to  //.  gives 


dd  dl3  da 

- b  2—  = - 

dji  dji  dji 


d<) 

dji 


did 


(h{ 


^  +  ^)  +  2—08  +  2—18 


djt 


dji 


da 


(4.73) 


At  the  critical  parameter  value  //  =  ji,c,  S  =  0,  j3  =  0,  and  7 2  =  c2  —  a.  Thus,  at 


d  dc- 


dS  dl3  da 

- b  2—  = - 

dji.  dji  dji 


dS 

dji. 


a) 


da 


Solving  these  equations  for  gives 


<13  1  c  da 

dji  2  c2  —  a  dji 


(4.74) 


(4.75) 


which  is  nonzero  if  ^7  7^  0  and  a  <  0. 

As  was  the  case  with  proposition  4.1,  the  final  step  in  the  proof  consists 
of  showing  that  all  other  eigenvalues  of  the  matrix  D  are  in  the  open  left  half 
complex  plane.  There  are  three  eigenvalues  of  D  which  correspond  to  one  nega¬ 
tive  real  value  eigenvalue  of  A  and  these  eigenvalues  are  solutions  of  the  (4.70). 
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By  using  the  Routh-Hurwitz  criterion,  we  can  show  that  solutions  of  equation 
(4.70)  in  C-  if  the  corresponding  real  eigenvalue  of  A  is  in  C'_ .  For  the  complex 
conjugate  pair  of  eigenvalues  of  A  (7,7),  we  have  following  two  equations 

A3  -  7A2  +  (c2  -  a) A  +  «7  =  0  (4.76) 

A3  —  7A2  +  (c2  —  a)A  +  ay  =  0  (4.77) 

Multiply  (4.76)  and  (4.77)  to  get  the  sixth  order  equation 

A6  -  (7  +  7)A5  +  (2(c2  -  a)  +  77)A4  +  (7  +  y)(2fl  -  c2)A3 

+  ((c2  —  a)2  —  2ft77)A2  +  (c2  —  «)«( 7  +  7)  A  +  ft2  77  =  0  (4.78) 

By  applying  the  Routh-Hurwitz  criterion  to  (4.78),  we  can  show  that  all  solu¬ 
tions  of  (4.78)  are  in  the  open  left  half  complex  plane  if  Re('y)  is  negative  (details 
are  in  Appendix  A).  □ 

We  have  proved  that  the  new  augmented  system  (4.65)  with  nominal  equi¬ 
librium  not  necessarily  at  the  origin  replaces  a  stationary  bifurcation  with  a 
codimension  two  bifurcation.  Note  that  the  system  has  the  same  critical  pa¬ 
rameter  value  for  the  original  system  and  the  augmented  system.  In  addition, 
the  crossing  eigenvalues  at  critical  point  are  located  0  and  ±\/c2  —  aj-  Also, 
note  that  an  original  simple  zero  eigenvalue  persists  under  the  augmentation. 
Thus,  monitoring  system  helps  by  introducing  juj  axis  eigenvalue  in  addition 
to  the  zero  eigenvalue.  Hence,  we  can  expect  that  the  power  spectrum  peaks 
near  bifurcation  to  be  located  at  0  and  V c2  —  a.  By  varying  c  and  a  (both 
tunable  parameters),  we  could  relocate  the  peak  at  \Jc2  —  a  to  any  desired  loca¬ 
tion.  This  flexibility  gives  more  assurance  to  say  that  peak  of  power  spectrum  is 
caused  by  closeness  to  instability  rather  than  other  factors  such  as  certain  noise 
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burst  at  that  specific  frequency.  However,  this  new  augmented  system  (4.65) 
also  comes  with  some  disadvantage  compared  to  the  system  Proposition  4.1.  In 
Proposition  4.1,  we  transform  a  stationary  bifurcation  into  a  Hopf  bifurcation. 
In  other  words,  the  system  has  a  periodic  orbit  as  its  solution  instead  of  a  new 
equilibrium  point  near  of  bifurcation.  In  comparison  to  the  previous  augmented 
system  design  (4.2),  the  new  augmented  system  (4.65)  shows  more  complicated 
bifurcation  behavior  [31].  The  system  is  no  longer  guaranteed  to  have  a  periodic- 
orbit  as  a  solution  near  bifurcation.  Either  a  periodic  orbit  or  a  new  equilibrium 
point  could  result  at  bifurcation.  The  bifurcation  diagram  depends  strongly  on 
the  vector  field  f(x.  //).  However,  it  may  be  possible  that  augmented  system  has 
desired  bifurcation  diagram  by  introducing  some  nonlinear  terms  into  augmented 
states.  Of  course,  to  do  that  we  have  detail  knowledge  on  f(x,p).  Details  on 
codimension  two  bifurcations  can  be  found  in  Guckenheimer  [31].  However,  for 
the  purpose  of  monitoring,  it  is  enough  to  have  a  discernible  power  spectrum 
peak  when  the  system  approaches  instability. 

The  next  proposition  asserts  that  the  new  augmented  system  also  works 
for  singular  perturbation  case  with  fewer  states  needed  to  be  fed  back.  Only 
difference  from  previous  singular  perturbation  case  is  that  we  no  longer  requires 
(HI)  of  Section  4.3.  Using  the  same  notation  in  the  Section  4.3,  we  have  the 
following  proposition. 

Proposition  4.7  Let  (H2)-(H5)  in  Section  f.3  hold  for  the  system  (f.l).  Then 
there  is  an  eo  >  0  and  for  each  eo  €  [0,  eo]  the  following  extended  system  undergoes 
a  codimension  2  (one  real  and  a  pair  of  complex  eigenvalues  crossing)  bifurcation 
at  an  equilibrium  m[[c’£  near  rn0  for  a  critical  parameter  value  qhc  near  /i.c. 

Xi  =  fi(x,  z,  fi,e)  —  cyi 
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iji  =  cxi  +  awi 

Wi  =  Vi 

ei  =  g(x.  z,  //.  e)  (4.79) 

where  i  =  1,  2, ....  n. 

Proof:  By  direct  application  of  Proposition  4.6  and  Theorem  7  of  [1].  □ 

We  are  going  to  just  state  following  proposition  without  detail  proof  for  the 
nonlinear  affine  control  system  since  it  is  direct  a  result  of  proposition  4.4  and 
proposition  4.6.  Here,  we  no  longer  need  assumption  (Bl)  which  assumes  origin 
remains  a  equilibrium  point  for  all  parameter  value  //,.  The  next  proposition 
refers  to  the  following  augmented  system: 

n 

X  =  f{x,(l)  +  J2di(X)Ui 

2=1 

y  =  c  JO  x  +  dz 

z  =  y  (4.80) 

where  x  €  Rn,  y  €  Rl,  c  €  i?|  and  /  is  l  x  l  identity  matrix.  We  have  new 
assumption  (B5)’  to  replace  assumption  (B5)  in  Section  4.3. 

(B5)’  The  subspace  D(x)  =  span{g.j{x)  :  i  G  l,...,m}  has  constant  dimension 
k  >  l  for  all:?:  in  Rn. 

Proposition  4.8  Assumptions  (B2)-(B4)  in  Section  4-3  and  (B5)’  hold  for 
the  system  (4-1),  then  the  augmented  system  (4-30)  experiences  a  codimension 
2  bifurcation  at  p  =  pc  with  special  input  u* .  In  addition,  if  the  nominal  equilib¬ 
rium  of  the  original  system  is  asymptotically  stable,  then  it  is  also  asymptotically 
stable  for  the  augmented  system  (4-80). 
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Proof:  Let’s  set  u*  =  —ch(x)y  where  h(x)  G  Rmxl  and  c  €  R.  Let  G(x)  = 
[gi{x),  ■  ■  ■ ,  gm.(x)\.  For  any  given  x,  G(x)  can  be  changed  into  singular  value 
decomposition  form  as  follows. 

G{x)  =  U{x)  X{x)V{x) 

(7l{x) 

0 
0 

0 

where  U(x)  and  V(x)  are  unitary  matrix,  i.e.,  Un (x)U(x)  =  I  and  VH (x)V(x)  = 
/,  with  dimension  n  x  n  and  m  x  m,  respectively.  Here,  A11  denotes  complex 
conjugate  transpose  of  A.  Because  of  assumption  (B5)’  and  smoothness  of  fj, , 
U(x),  V(x),  and  E(x)  are  smooth  function  of  a:.  Set  h(x)  to 


u\(x) 


un{x) 


0 

a2(x) 

0 


0 

0  0 
Vl 

0 


v\{x) 

Vm{x) 


(4.81) 


h{x)  =  VH{x  )M(x) 


E(x)M(x) 


uf  ■  ■  ■  uj 1  0  •  •  •  0 


(4.82) 


Then, 


G{x)h{x) 


I 

0 


(4.83) 


where  /  is  l  x  l  identity  matrix.  Rest  of  the  proof  is  direct  application  of  Propo¬ 
sition  4.4.  □ 
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Chapter  5 


Application  of  Monitoring  System  to  Axial 
Flow  Compression  System 


There  have  recently  been  several  important  developments  in  analysis  and  control 
of  axial  flow  compressor  both  in  analysis  of  instability  phenomena  and  their  con¬ 
trol.  These  developments  make  possible  the  increased  performance  of  the  axial 
flow  compressor  by  operating  near  the  maximum  pressure  rise.  However,  the 
increased  performance  comes  at  the  price  of  significantly  reduced  stability  mar¬ 
gin.  Loss  of  stability  is  of  course  unacceptable.  It  results  decreased  performance 
of  the  compressor  and  to  mechanical  damage  of  the  compression  system.  Axial 
flow  compressor  are  known  to  be  susceptible  to  two  basic  types  of  instability, 
[29].  One  of  these  is  surge  which  is  a  low-frequency,  large-amplitude  oscillation 
of  the  mean  mass  flow  rate  and  pressure  rise.  The  other  is  rotating  stall  which 
corresponds  to  a  traveling  wave  of  gas  around  the  annulus  of  the  compressor. 
These  rich  nonlinear  characteristics  of  the  axial  flow  compressor  are  well  suited 
to  apply  our  suggested  monitoring  system. 

Several  control  laws  have  been  introduced  to  permit  operation  near  peak 
pressure  rise  while  maintaining  stability  (e.g.  [35],  [40],  [51],  [56],  etc.).  Many  of 
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those  suggested  control  laws  ([50]  and  [51])  are  employ  linear  control  for  avoiding 
or  delaying  the  occurrence  of  stall.  There  are  also  nonlinear  controls  ([11],  [56], 
[39]  and  [40])  that  aims  to  ensure  results  in  only  stable  rotating  stall.  Thus,  even 
though  the  nominal  equilibrium  is  not  stabilized,  it  may  be  possible  to  stabilize 
a  neighborhood  of  the  nominal  solution  for  a  range  of  parameter  values  includ¬ 
ing  the  stall  value  of  the  throttle  opening  parameter,  to  finite  perturbations. 
However,  there  is  no  known  control  law  which  totally  eliminates  the  bifurcation 
and  maintains  stable  nominal  equilibrium  point  for  all  parameter  values.  Since 
bifurcation  will  occur  for  any  controller,  the  axial  flow  compressor  is  a  good 
application  for  our  monitoring  systems. 

5.1  Modeling  and  Bifurcation  Analysis 

Compressor  dynamic  modeling  is  a  subject  that  has  attracted  significant  atten¬ 
tion  (see,  e.g.,  [20],  [24],  [14]  and  [33]).  Greitzer  [30]  developed  a  nondimensional 
fourth  order  compression  system  model  and  introduced  a  nondimensional  pa¬ 
rameter,  B ,  which  he  found  to  be  a  determinant  of  the  nature  of  post-instability 
behavior.  A  global  bifurcation  of  periodic  solutions  and  other  bifurcations  were 
found  for  this  model  [41],  and  were  used  to  explain  the  observed  dependence  of 
the  dynamical  behavior  on  the  B  parameter. 

Moore  and  Greitzer  [48]  introduced  a  refined  model  to  describe  stall  phenom¬ 
ena  in  axial  compressor.  This  model  includes  the  dynamics  of  nonaxisymmetric 
flow  patterns,  which  was  not  present  in  the  Greitzer  model  [30]. 

For  this  dissertation,  we  employ  a  Moore  and  Greitzer’s  model  in  [48]  with 
slight  modification.  Moore  and  Greitzer  derived  following  basic  partial  differen- 
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Figure  5.1:  Schematic  of  an  axial  compressor 

tial  equation  representing  a  local  momentum  balance  in  the  compressor  and  its 
associated  ducting: 

Ap  =  C„(V  +  vo)  -  lc%  -  <4  4  vdn  -  (5.1) 

where  V  denotes  the  annulus-averaged  (mean)  gas  axial  velocity;  Vo  is  the  axial 
velocity  perturbation  evaluated  at  r/  =  0  (the  inlet  face  of  the  compressor);  AP 
is  the  plenum  to  atmosphere  pressure  rise;  r ),  6  are  the  axial  and  angular  coordi¬ 
nates,  respectively;  and  a ,  lc.  m  are  internal  compressor  lag,  overall  compressor 
length,  and  exit  duct  length  factor,  respectively. 

The  compressor  characteristic  Css  is  particular  to  each  compressor.  Moore 
and  Greitzer  [48]  used  following  cubic  function  as  a  compressor  characteristic 

C«(  14*)  =  /0  +  1/[1  +  |(—  -  1)  -  h—  -  l)3]  (5.2) 

Z  UJ  Z  UJ 

where  \]oc  =  V  +  vq  (the  total  local  axial  flow);  fo  is  shut-off  head;  w  is  a 
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compressor  characteristic  width  factor;  and  H  is  a  is  a  compressor  characteristic 
height  factor.  Unlike  Moore  and  Greitzer’s  model,  we  only  assume  here  that  Css 
be  a  smooth  function  of  V  +  v0  as  in  Liaw  and  Abed  [40].  If  there  are  no  spatial 
variation  of  gas  density  and  pressure  in  the  plenum,  an  overall  material  balance 
on  the  gas  over  plenum  gives: 

^f  =  1L[V-F(7,AP)J  ,5.3) 

where  F(y,  AP)  is  a  inverse  of  throttle  characteristic  and  7  is  a  parameter  pro¬ 
portional  to  the  throttle  opening.  Note  that  parameter  B  which  was  first  intro¬ 
duced  by  Greitzer  [30]  appears  in  (5.3). 

Moreover,  Moore  and  Greitzer  retain  only  the  first  harmonic  of  the  ax¬ 
ial  velocity  perturbation  and  assume  that  a  first  harmonic  axial  velocity  per¬ 
turbation  is  a  sine  wave  with  time  varying  phase  and  amplitude  (i.e.  «o  = 
W  A(t)sin{6  +p{t))).  By  employing  Galerkin  procedure  to  one  mode  truncation 
approach,  the  model  is  reduced  to  the  following  3  dimension  nonlinear  ODE: 

dA  n  f 271 

—  =  CSS(V  +  W A  sin  0)  sin  6cW  (5.4) 

clt  ttW  Jo 

dV  1 

—  =  -A  P  +  —  CSS(V  +  W  A  sin  0)d0  (5.5) 

dt  27 r  Jo 

=  ^lv-F(r^P)]  (5.6) 

Stability  analysis  has  been  done  for  compressor  without  employing  bifurca¬ 
tion  theory  (  [19], [44]  and  [49]).  McCauglian  [43]  performs  bifurcation  analysis 
on  model  of  Moore  and  Greitzer  [48]  and  shows  that  a  stationary  bifurcation  from 
nominal  equilibrium  point  occurs  as  the  throttle  opening  parameter  7  is  varied. 
Liaw  and  Abed  [40]  extended  the  work  of  McCauglian  [43]  to  the  case  in  which 
the  axisymmetric  compressor  characteristic  is  taken  to  be  a  general  smooth  func- 
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tion  of  axial  velocity.  Since  Liaw  and  Abed  [40]  did  bifurcation  analysis  of  the 
axial  compressor  model  (5.4)-(5.6),  we  follow  the  bifurcation  analysis  therein. 

We  assume  that  B  >  0  and  F  is  strictly  increasing  function  with  respect  to 
each  of  the  variables  7  and  A P  in  (5.6).  Seeking  equilibrium  point  of  (5.4)-(5.6), 
we  note  that  A  =  0  always  results  in  =  0.  There  may  also  be  equilibrium 
points  for  which  A  ^  0.  Consider  a  nominal  equilibrium  point  for  which  there 
is  no  asymmetry  flow  ,  i.e.  let  A  =  0.  Denote  this  nominal  equilibrium  point  by 

r  1 T 

x°h)  =  0  C°(7 )  A  P°(7)  (5.7) 

where  and  A P°(7)  satisfy  V°  =  F(7,  AP°)  and  A P°  =  CSS{V°).  At  this 

equilibrium  point,  the  system  Jacobian  matrix  is 

aV„(r°(7))  0  0 

£0=  0  C'„(l’“(7))  -1  (5-8) 

_0  -^DzrF(7,AP°)  _ 

At  the  parameter  value  7  =  7°  for  which  Cfaa{VQ{^))  =  0,  L0  has  a  zero  eigenvalue 
and  two  eigenvalue  with  negative  real  part  since  we  assume  that  B  >  0  and  F 
is  strictly  increasing  function  with  respect  to  each  of  variables  7  and  AP.  This 
suggests  that  a  static  bifurcation  occurs  at  7  =  70.  Moreover,  it  is  not  difficult 
to  see  that  for  the  matrix  L0  to  have  a  pair  of  pure  imaginary  eigenvalues,  C'ss 
must  be  positive.  Because  the  trace  of  the  lower  right  2x2  submatrix  of  L0 
should  be  zero.  However,  if  this  were  the  case  then  the  matrix  /,0  would  have 
aC'ss(V°h))  as  a  positive  real  eigenvalue,  and  the  equilibrium  would  therefore 
unstable.  This  means  that  a  stationary  bifurcation  must  be  occur  before  an  Hopf 
bifurcation  for  the  nominal  equilibrium  of  ,t°.  In  the  remainder  of  this  section, 
we  focus  on  the  stationary  bifurcation  that  occurs  for  7  =  70. 
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To  analyze  bifurcation  behavior  of  the  model,  we  employ  bifurcation  formulae 
in  Section  2.3.  Let  .t0  be  equilibrium  point  at  which  C^S(V0(7))  =  0  for  7  =  70. 
The  Taylor  series  expansion  of  (5.4)-(5.6)  for  (x,j)  near  (x°,7°)  is  given  by 


dx 

—  =  L„x  +  Qo(x,x)  +  C„(x,x,x)  +  (7  -  i  )Ltx 


(5,9) 


where  x  now  denotes  the  deviation 
and 


1  T 


A  V  Ar 


—  ,t°.  Here,  L0  is  as  in  (5.8), 


/ 


Qo(x,x)  = 


aC”s(V°)x \x2 

lC:s(V°)(W2xl  +  2xl) 


(5.10) 


V  -^%r)^(70,AP°)i1  ) 

\aC"s's{Va){w2x\  +  Axixl) 


(  1 


Cq(x,x,x)  = 


25  \ 


C^V^W^x  2  +  Ixl) 


(5.11) 


Li  = 


\ 


V  -^D(ap)sF(^,AP°)xI  ) 

O'0)6  0  0 

0  C;s(H°)6  0  (5.12) 

0  0 

where  £1  =  D1F{rf,  AP°)  and  £2  =  -Dap7-F(7°,  AP°).  Set  1  =  10  0  and 

r  =  lT\  the  left  and  right  eigenvectors  respectively  corresponding  to  the  zero 
eigenvalues  of  L0.  Check  the  transversality  condition  by  using  (2.43) 


( 


V 


lL1r  =  aC:s(V°)D1F(1\AP°) 


(5.13) 


We  calculate  the  bifurcation  stability  coefficients  using  (2.44)  and  (2.45): 


Pi  =  lQo(r,r)=0  (5.14) 

p2  =  2/{2<20(r,  x2)  +  C0(r,  r,  r)} 

=  ±aW2{2DAPF(%  AP°)[C£(V“)]2  +  CftV*)}  (5.15) 
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The  next  theorem  follows  from  the  foregoing  discussion. 


Theorem  5.1  Suppose  that  C"ss  ^  0  and  that  F  is  strictly  increasing  in  each 
of  its  variables.  Then  the  axial  compression  system  model  (5.f)-(5.6)  exhibits  a 
pitchfork  bifurcation  with  respect  to  small  variation  of  7  at  the  point  (a’°,7°)  for 
which  CSS{V°{^))  =  0. 

Note  that  the  formula  (5.15)  shows  that  stability  of  bifurcation  equilibria 
depends  on  derivatives  of  the  axisymmetric  compressor  characteristic  at  the  bi¬ 
furcation  point. 

To  illustrate  Theorem  5.1  numerically,  we  consider  an  example  in  which  the 
compressor  characteristic  is  cubic.  Let  the  axisymmetric  compressor  character¬ 
istic  CSS(V)  and  the  throttle  characteristic  F  be 

CSS{V)  =  1.56  +  1.5(V  —  1)  -0.5(F-  l)3  (5.16) 

F(y,  AP)  =  7VA P  (5.17) 

Note  that  (5.16)  is  the  same  cubic  function  which  used  in  Moore  and  Greitzer 
[48]  (see  (5.2)).  Also,  throttle  characteristic  (5.17)  is  same  throttle  character¬ 
istic  given  in  Moore  and  Greitzer  [48]  and  satisfies  strictly  increasing  function 
condition. 

Select  parameter  values  a  =  3 ,W  =  1.0,  and  B  =  0.5.  By  using  the  bifur¬ 
cation  analysis  package  AUTO  [22],  we  generate  Figure  5.2.  Figure  5.2  shows 
the  first  bifurcation  occurring  at  7  =  1.25  and  second  bifurcation  occurring  at 
7  =  1.17325.  First  one  is  a  pitchfork  bifurcation  corresponding  to  rotating  stall 
and  second  one  is  Hopf  bifurcation.  In  the  Figure  5.2,  solid  line,  dotted  line,  and 
circle  correspond  to  stable  equilibrium  point,  unstable  equilibrium  point,  and 
periodic  orbit,  respectively.  Figure  5.2  (b)  is  the  blown  up  figure  of  Figure  5.2 
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(a)  for  7  value  between  and  1.3  where  fill  up  circle  corresponds  to  stable  periodic- 
solution  and  empty  circle  corresponds  to  unstable  periodic  solution.  This  figures 
are  drawn  for  variable  A. 


5.2  Precursor  for  Rotating  Stall 

In  the  previous  section,  we  have  shown  that  the  axial  compression  system  model 
(5.4)-(5.6)  experiences  a  pitchfork  bifurcation  as  the  parameter  7  varies.  Based 
on  this  information,  we  will  apply  our  monitoring  system  (4.2)  in  Section  4.1  to 
the  axial  flow  compressor. 

The  Jacobian  matrix  of  the  compressor  model  (5.8)  shows  that  the  linearized 
system  is  decoupled  as  assumption  (B3)  of  Proposition  4.4  between  A  and  V, 
A P.  Thus,  we  need  only  one  augmented  state  corresponding  to  A.  However, 
to  apply  our  monitoring  system,  we  should  have  monitoring  input  to  state  A. 
But,  our  compressor  model  ((5.4)-(5.6))  does  not  have  input  to  amplitude  of 
first  harmonic  of  asymmetric  flow  (A).  Therefore,  we  are  going  to  adopt  new 
compressor  model  which  has  control  input  to  state  A. 

It  is  known  that  there  are  many  ways  to  generate  asymmetric  flow  in  an  axial 
compressor,  such  as  oscillating  the  inlet  guide  vanes,  vanes  with  oscillating  flaps, 
jet  flaps,  tip  bleed  above  the  rotor,  etc.  In  Paduano  et  al.  [51],  control  was 
implemented  using  a  circumferential  array  of  hot  wires  to  sensor  propagating 
waves  of  axial  velocity  upstream  of  the  compressor.  Using  this  information, 
additional  circumferential  traveling  waves  were  then  generated  with  appropriate 
phase  and  amplitude  by  wiggling  inlet  guide  vanes  driven  by  individual  actuator. 
Their  modified  compressor  schematics  is  shown  in  Fig  5.3. 
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Figure  5.3:  Schematic  of  controlled  IGV  axial  compressor 


In  another  paper  [50],  Paduano  et  al.  obtained  a  linear  state  space  model 
which  includes  the  effect  of  the  input  (moving  the  inlet  guide  vanes).  In  [50], 
experiments  were  used  to  identify  the  effect  of  input.  Until  now,  no  low  order 
nonlinear  model  is  available  that  includes  the  effect  of  the  dithering  the  inlet 
guide  vanes.  Therefore,  we  settle  for  using  the  following  general  nonlinear  model 
in  this  section.  Let  the  compression  system  with  controlled  inlet  guide  vanes  be 
described  by 


clA 

dt 

dV 

dt 

dAP 

dt 


a 

7lTF 


r  2-7T 


CSS{V  +  W A  sin  6)  sin  6d6  +  gi{A,  V,  A P)u 


1  r2n 

-A P  +  —  /  Caa(V  +  W A  sin  d)dd  +  g2{A ,  V,  A P)u 
2i r  Jo 


4  B2 


[in,-  F(~.AP)[  •  g,(A.\.  Al>)u 


(5.18) 


Here,  we  assume  that  the  effect  inputs  to  the  system  can  be  modeled  as  affine. 
However,  since  dynamics  of  AP  is  derived  from  balance  of  entering,  leaving,  and 
stored  mass  of  the  plenum,  changing  the  phase  and  amplitude  of  traveling  wave 
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A  is  most  likely  not  affect  the  dynamics  of  AP. 

We  have  the  following  proposition  for  the  model  above. 


Proposition  5.1  Assume  that  gi,g2,  and  g 3  are  general  smooth  function  with 
respect  to  each  of  their  variables  and  gi  does  not  depend  on  the  bifurcation  pa¬ 
rameter.  In  addition,  suppose  that  gi  does  not  change  sign  along  the  equilibrium 


path  as  the  parameter  7  varies  and  its  sign  is  k7iow7i  to  us  but  not  the  exact  func¬ 
tion  itself.  Then  following  augmented  system  with  u  =  —sgn(gi)cy  undergoes  a 
Hopf  bifurcation  at  70  which  corresponds  to  a  pitchfork  bifurcation  point  from 


nominal  equilibrium  point  x°  (5.7)  for  the  original  system: 


(I  A  ry  r2ir 

—  =  —  Css{  V  +  WA  sin  9 )  sin  9d6  +  gA A,  V,  A P)u 

dt  7T  14  J  0 

=  cA 
dt 

dV  1  r2* 

—  =  -AP  +  —  /  CSS{V  +  WA  sin  9)d9  +  g2{A ,  V,  A P)u 

dt  Ztt  Jo 


dAP 

dt 


:[V-F(j,AP)]  +  g3(A.V,AP)u 


(5,19) 


Proof:  Linearize  the  equation  (5.19)  at  the  equilibrium  point.  The  equilibrium 
point  for  system  (5.19)  corresponding  to  a  nominal  equilibrium  point  x°  is  given 
by 

r  1 T 

0  0  P°(7)  AP°(7)  (5.20) 

Note  that  augmentation  does  not  change  the  nominal  equilibrium  point  x°.  Then 
the  Jacobian  matrix 

aCjV°(7))  -s  |  Si  |  0  0 

c  0  0  0 

(5.21) 

0  92  C„(V»(7))  -1 

0  <]3  —j^D^F^jAP0) 


79 


where  r/2  :=  g2sgn(g\)  and  g?,  :=  gzsgn(gi).  Since  the  Jacobian  matrix  is  in 
block  lower  triangular  form,  we  can  apply  Proposition  4.1  to  the  upper  left  2x2 
submatrix  to  verify  occurrence  of  Hopf  bifurcation.  At  the  critical  point,  we 
have  a  pair  of  imaginary  eigenvalues 


A 


±c\/l  gi(')  1 3 


(5.22) 


□ 


From  equation  (5.22),  without  knowing  g\  we  cannot  determine  exact  location 
of  crossing  pair  of  imaginary  eigenvalues.  However,  this  equation  does  show  that 
the  values  of  the  eigenvalues  at  crossing  are  proportional  to  the  value  of  c,  which 
we  can  control. 

Now  consider  the  case  in  which  g\  is  exactly  known.  We  have  the  following 
Proposition. 

Proposition  5.2  Consider  the  augmented  system  (5.19).  Assume  that  g j-1  ex¬ 
ists  and  is  independent  of  the  bifurcation  parameter  7.  Then  the  system  (5.19) 
with  u  =  —gilcy  undergoes  a  Hopf  bifurcation  for  7  =  70  which  corresponds 
to  a  pitchfork  bifurcation  point  for  the  original  system.  Also,  a  crossing  pair  of 
imaginary  eigenvalues  has  value  ±ci  at  the  critical  point. 

Note  that  if  gf1  fails  to  exist  in  the  region  of  interest,  then  we  can  use  Proposition 
5.1.  However,  we  can  determine  the  exact  location  of  crossing  pair  of  imaginary 
eigenvalues  unlike  Proposition  5.1  since  we  have  exact  knowledge  of  g\. 

In  the  rest  of  this  section,  we  try  to  verify  the  propositions  above  numeri¬ 
cally.  All  of  the  bifurcation  diagram  shown  here  are  obtained  using  numerical 
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bifurcation  analysis  package  AUTO.  We  assume  that  g\  and  g2  are  given  as 


gi(V)  =  v2  + 1 
g2(AP)  =  A  P2 

93  =  0 

We  also  take  the  same  compressor  characteristic  (5.16),  inverse  function  of  throt¬ 
tle  pressure  rise  (5.17),  and  parameters  which  are  given  in  Section  5.1.  Moreover, 
small  white  Gaussian  noise  (ni(£))  is  added  to  ^  in  equation  (5.18)  either  nat¬ 
urally  or  artificially  such  that 
r]  A  n  A77 

—  =  Cas(V  +  WA  sin  6)  sin  OdO  +  gx (A,  V,  A P)u  +  m{t)  (5.23) 

dt  n  W  Jo 

This  assumes  that  there  exist  noise  disturbance  in  the  natural  environment.  We 
have  shown  that  the  compressor  with  these  specifications  undergoes  pitchfork 
bifurcation  at  7  =  1.25  without  input. 

Consider  problem  setting  of  Proposition  5.2.  Since  g\  is  known,  we  can  cancel 
out  (]\  in  (5.23)  and  result  g\U  =  —  cy  by  setting  input  u  =  —  v21+1  cy.  Figure  5.4 
(a)  shows  that  augmented  system  undergoes  Hopf  bifurcation  and  the  bifurcation 
point  is  the  same  as  original  system  parameter  value.  In  Figure  5.5  where  we 
set  c  =  5,  note  that  the  clear  pronounced  spectrum  peak  is  at  5  as  7  approaches 
critical  value.  Another  Figure  5.6  shows  that  the  power  spectrum  peak  is  moved 
to  10  when  c  is  changed  to  10. 

Let’s  consider  the  case  such  that  there  is  no  way  to  add  artificial  disturbance 
to  A  and  also  no  natural  disturbance  in  A.  Since  existence  of  a  small  white 
Gaussian  noise  is  crucial  to  see  a  growing  peak  in  the  power  spectrum  as  a 
system  closes  to  a  bifurcation,  we  have  to  add  some  noise  to  the  system.  Here, 
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Power  Spectral  Density 


x  1  O  4 


a.  7  =  10 


5  lO  15  20  25  30 

Frequency  (rads/sec) 


b.  7=  1.3 


Figure  5.5:  Power  spectrum  of  A  when  c=5 


we  are  going  to  add  a  white  Gaussian  noise  {ri2(t))  to  the  augmented  state  y 
such  that 

^  =  CA  +  n2(t)  (5.24) 

Since  we  have  total  control  over  the  augmented  state  y,  it  is  easy  to  add  noise 
to  state  variable  y. 

With  noise  injected  into  the  y  dynamics  as  in  (5.24),  Figure  5.7  shows  that 
we  have  a  peak  at  co  =  5  for  power  spectrum  measurement  of  both  A  and  y. 
This  shows  that  even  if  noise  does  not  enter  naturally  or  cannot  be  injected  into 
the  physical  system  dynamics,  it  may  still  be  possible  to  obtain  a  precursor  by 
injecting  noise  into  the  augmented  states. 

Next,  we  consider  Proposition  5.1  case.  This  case  we  cannot  cancel  out  gi 
since  we  assume  no  knowledge  on  gL  except  its  sign  which  is  locally  non  varying. 
Figure  5.4  (b)  shows  that  augmented  system  still  undergoes  Hopf  bifurcation 
at  the  same  critical  value  despite  of  partially  known  g\  assumption.  Let  our 
monitoring  input  u  =  —  cy  where  c  >  0  because  gi(V)  =  V2  +  1  >  0  for  all  mc 
In  the  Figure  5.8  where  we  set  c  =  5,  we  see  clear  power  spectrum  peak  around 
7  which  we  could  predict  from  (5.22)  as  7  closes  its  bifurcation  point.  Note  that 
power  spectrum  peak  moves  to  around  14.14  by  changing  c  to  10  in  the  Figure 
5.9. 

5.3  Precursor  for  Hopf  Bifurcation  of  the 
Axisymmetric  Solution 

From  previous  sections,  we  know  that  axial  compressor  undergoes  a  pitchfork 
bifurcation  as  parameter  varies.  Here,  we  are  going  to  consider  the  case  of  after 
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Frequency  (rads/sec) 


a.  7  =  1.3  and  power  spectrum  of  A 


5  1  O  15  20  25 

Frequency  (rads/sec) 


b.  7=  1.3  and  power  spectrum  of  y 


Figure  5.7:  Power  spectrum  with  noise  driving  only  the  augmented  state  y 


pitchfork  bifurcation  occurred.  For  this  section,  we  limit  our  case  to  a  cubic 
compressor  characteristic  which  is  given  by  equation  (5.16).  Also,  we  employ 
throttle  characteristic  given  by  equation  (5.17). 

After  the  pitchfork  bifurcation  occurs,  there  are  two  new  equilibrium  points 
for  the  amplitude  of  the  first  harmonic  of  asymmetric  flow  (A).  Moreover,  the 
equilibrium  point  with  A  =  0  is  no  longer  stable.  However,  if  initially  A  =  0  and 
there  is  no  perturbation  to  dynamics  of  A ,  then  A  remains  zero  for  all  the  time 
(see  equation  (5.4)).  This  corresponds  to  an  axisymmetric  flow  condition  and 
it  is  stable  for  given  cubic  compressor  characteristic  and  throttle  characteristic. 
For  this  ideal  compressor  which  A  being  fixed  to  0  even  after  the  pitchfork 
bifurcation,  the  first  bifurcation  will  be  a  Hopf  bifurcation. 

Since  we  assume  that  A  remains  zero  all  the  time,  we  can  drop  equation  (5.4) 
and  the  following  two  dimensional  ODE  remains: 


dV_ 

dt 

dAP 

dt. 


-A  P  +  C8S(V) 


1 

~\1P 


[V-F(7,AP)] 


(5.25) 


Figure  5.2  shows  the  Hopf  bifurcation  occurring  at  7  =  1.17325  for  system  (5.25). 

Before  we  apply  our  monitoring  system,  there  are  two  things  to  note.  First, 
the  equilibrium  point  varies  as  the  parameter  7  varies.  Hence,  we  are  going 
to  use  our  monitoring  system  with  wash  out  filter  (4.65).  Second,  we  have  to 
have  monitoring  signal  input  to  both  states  of  the  system  (5.25).  However,  our 
axial  compressor  has  one  input  as  form  of  7  =  7"  +  u.  From  this,  it  is  clear 
that  we  cannot  transform  axial  compressor  model  system  directly  into  suggested 
augmented  state  form.  Therefore,  we  suggest  a  modified  method  as  follows. 
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The  modified  augmented  model  is 


dV_ 

dt 

dAP 

dt 

dy 

dt 

dz 

dt 


-A  p+±CS8(V) 

Z7T 

jfjz.v  ~  7n^A P]  -  cy 
cAP  +  dz 

y 


(5.26) 


We  set  7  =  7"  +  1v/~7p  cy  to  get  above  equations  where  7"  denotes  original 
parameter  setting.  Note  that  equilibrium  point  of  (5.25)  is  still  the  equilibrium 
point  of  (5.26). 

The  linearization  of  (5.26)  at  equilibrium  point  ( V'°,  AP’O,  — |A P°)  where  V° 
and  A P°  satisfy  V°  =  P(y,  AP°)  and  A P°  =  6'.,., ( T  '° ) .  gives  Jacobian  matrix 


D  = 


C'JV°) 

1 

IB2 

0 

0 


-1 

ryn 


0  0 
-c  0 


8BVAP0 

c  0  d 

0  1  0 


The  characteristic  polynomial  of  D  is 


(5.27) 


A 1  —  (a  +  (3) A3  +  (ck/J  +  q(g2  T  —  /) A2  +  ( (ck  +  (3)d  —  ckc2) A  —  d(c k/3+  )  (5.28) 


where  tk  =  C'ss  ( T'° )  and  (3  =  —  8b2  To  check  the  stability  of  augmented 

system  if  all  the  eigenvalues  of  Jacobian  matrix  of  the  system  (5.25)  is  in  the 
open  left  half  complex  plane,  we  to  use  the  Routh-Hurwitz  criterion.  The  Routli 
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array  is 


?3  —(a  +  {3) 

2  (“+/3)(afl+-j-p')+,3c2 

a+f. 3 

I  c2(-Q(afl+-j^2-)(Q+,!3)-Qfl(c2-^)+/32rf) 
(a+/3)(a/3+^)+/3c2 

5°  -^0^+352) 


ci/J  +  4^5  'I  r2  (I  —cl(a/3  +  yp) 

((a  +  /3)d  —  ac2)  0 

—d(a[3  +  0 

0  0 


0 


0 

(5.29) 


Note  that  (a  +  (3)  <  0  since  original  system  has  eigenvalues  in  open  left  half 
complex  plane  and  aft + p?  >  0  since  determinant  of  Jacobian  matrix  of  (5.25)  is 
positive  if  (5.25)  has  all  the  eigenvalues  in  the  open  left  half  complex  plane.  From 
these,  we  have  following  sufficient  conditions  to  guarantee  asymptotic  stability 
of  augmented  system. 


ft  <  0  and  d  <  —  \K 


(5.30) 


where  d  is  constant  which  we  can  choose  and  should  be  chosen  big  enough  to 
guarantee  the  stability  (see  s-1  row).  Since  3  =  —  @ 's  a'wabs  less  than 

zero  if  the  equilibrium  point  of  (5.26)  is  real  and  7"  is  positive.  These  are  the 
case  for  our  axial  compressor  set-up. 

By  using  the  fact  that  at  criticality  a  +  (3  =  0,  the  characteristic  polynomial 
at  the  criticality  simplifies  to 


A4  +  (cj2  +  c2  -  d) X2  -  ac2 X  +  duj2  =  0 


(5.31) 


where  lo  =  ± 


p-  is  a  crossing  pair  of  imaginary  eigenvalues  of  (5.25). 


4B2 


For  the  compressor  with  cubic  characteristic  function,  a  is  close  to  zero  near  the 
critical  point.  Hence,  we  can  expect  equation  (5.31)  has  solutions  close  to 
\jd  —  (c2  +  uj2)  ±  yj d 2  —  2 dc  +  c4  +  2du)2  +  2 c2u)2  +  co4 


A  =  ±- 
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(5.32) 
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It  is  easy  to  see  that  inside  of  outer  square  root  of  above  equation  is  always 
less  than  zero.  Therefore,  augmented  system  has  two  pair  of  complex  conjugate 
eigenvalues  near  critical  point. 

As  we  have  pointed  out  before,  the  augmented  system  has  two  pair  of  complex 
conjugate  eigenvalues  with  real  part  guaranteed  to  be  less  than  zero  but  very 
close  to  imaginary  axis  near  bifurcation  point.  Actual  bifurcation  diagram  is 
not  shown  here  since  until  now  there  are  no  numerical  package  which  detects 
codimension  2  bifurcation.  However,  Figure  5.3  which  is  generated  by  AUTO 
shows  stability  change  in  equilibrium  point  at  same  critical  value  of  7  as  original 
system. 

Despite  of  this  unclear  post  bifurcation  behavior,  we  can  still  use  augmented 
system  to  get  more  assurance  for  closeness  to  bifurcation  point  as  follows.  First, 
we  observe  the  power  spectrum  without  augmentation.  Once  we  saw  power 
spectrum  peak  at  lo.  then  apply  feedback  with  state  augmentation.  If  the  power 
spectrum  peak  moves  the  point  where  it  supposed  to  be,  then  we  have  more 
assurance  that  system  is  really  close  to  its  instability. 

Note  that  original  system  has  a  complex  conjugate  pair  of  eigenvalues  ±i  at 
bifurcation  point  for  our  chosen  parameter.  From  the  calculation  of  section  3.2, 
we  can  expect  power  spectrum  peak  at  co  =  l  and  we  can  verify  that  from  Figure 
5.11. 

For  the  augmented  system,  we  see  two  peaks  as  system  close  to  its  bifurcation 
point.  In  Figure  5.12  where  we  have  chosen  c  =  5  and  cl  =  —100,  we  see  to  peaks 
around  at  lo  =  1  and  co  =  11.  We  have  moved  these  power  spectrum  peak 
to  lo  =  1  and  lo  =  15  by  varying  d  to  -200  in  Figure  5.13.  Those  location  of 
power  spectrum  peak  could  be  predicted  by  solving  (5.31).  We  can  say  that 
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our  augmentation  does  not  have  much  effect  on  original  power  spectrum  peak 
location  rather  we  create  new  power  spectrum  peak  which  is  more  influenced  by 
changing  cl  than  original  peak. 
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Chapter  6 


Conclusions  and  Suggestions  for  Future 
Research 


We  have  discussed  that  a  power  spectrum  precursor  can  be  a  useful  measurement 
to  detect  nearness  to  instability  in  a  system.  It  has  been  shown  that  a  power 
spectrum  precursor  is  a  useful  signal  of  closeness  to  bifurcation  in  the  case  of  Hopf 
bifurcation.  Unlike  Hopf  bifurcation,  employing  a  power  spectrum  precursor  to 
stationary  bifurcation  does  not  result  in  satisfactory  performance.  To  circumvent 
this  we  have  suggested  an  augmented  monitoring  system. 

The  augmented  system  has  been  shown  to  keep  most  of  the  original  nonlinear 
characteristics  of  the  system,  except  that  it  replaces  a  stationary  bifurcation  by 
a  Hopf  bifurcation.  The  transformed  system  undergoes  bifurcation  at  the  same 
critical  paremeter  value  as  the  original  system.  Moreover,  it  has  been  shown  that 
we  can  transform  either  supercritical  and  subcritical  pitchfork  bifurcation  to  a 
supercritical  Hopf  bifurcation  under  mild  assumptions.  This  achieves  both  pre¬ 
cursor  for  bifurcation  and  stabilizing  the  bifurcated  solution.  Moreover,  should 
the  system  undergo  a  Hopf  bifurcation,  the  augmented  system  undergoes  a  codi¬ 
mension  two  bifurcation.  Thus  the  approach  allows  detection  of  both  stationary 
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and  Hopf  bifurcations. 

One  advantage  of  the  augmented  system  is  that  it  does  not  require  detailed 
knowledge  of  the  system  being  monitored.  Mere  knowledge  that  the  system 
undergoes  stationary  bifurcation  is  enough  to  rendering  stationary  bifurcation 
to  Hopf  bifurcation.  Another  advantage  is  that  we  have  to  total  control  on 
the  frequency  of  bifurcated  periodic  solution.  This  means  we  can  only  measure 
certain  range  of  frequency  rather  than  measuring  the  entire  frequency  domain. 

These  advantages  are  achieved  under  some  restrictive  assumptions  such  as  the 
origin  being  an  equilibrium  point  as  the  parameter  varies  and  having  independent 
inputs  into  all  original  system  states.  We  have  discovered  that  in  some  cases  we 
can  transform  stationary  bifurcation  to  Hopf  bifurcation  with  restricted  input. 
Moreover,  for  singularly  perturbed  system  case  we  need  inputs  only  to  the  slow 
varying  states.  For  the  case  of  non  origin  equilibrium  point,  we  used  wash 
out  type  filter.  However,  it  renders  stationary  bifurcation  to  codimension  2 
bifurcation. 

We  have  successfully  applied  a  precursor  to  axial  compressor.  Numerical 
simulation  on  compressor  have  shown  that  we  have  a  clear  precursor  as  a  system 
closes  a  bifurcation.  We  have  examined  a  precursor  for  both  of  a  stationary 
bifurcation  corresponding  to  rotating  stall  and  a  Hopf  bifurcation  corresponding- 
surge.  It  should  be  noted  that  with  some  detail  knowledge  of  the  system  we 
can  do  more.  This  is  certainly  the  case  for  axial  flow  compressors,  especially  for 
surge  detection. 

Several  directions  for  extension  of  this  work  are  summarized  as  follows. 

1.  In  the  derivation  of  a  power  spectrum  precursor,  we  assumed  that  all  other 
eigenvalues  have  relatively  large  negative  real  part  compared  to  imaginary 
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axis  crossing  critical  eigenvalue.  Of  course  some  systems  might  have  two  or 
more  eigenvalues  near  the  imaginary  axis.  Moreover,  in  some  cases  system 
could  have  imaginary  axis  approaching  eigenvalues  without  ever  crossing 
imaginary  axis  as  parameter  varies.  A  power  spectrum  precursor  could  not 
discern  these  cases.  Research  to  address  these  issues  is  recommended. 

2.  We  assumed  that  white  Gaussian  noise  enters  the  dynamical  model  in  the 
form  of  additive  forcing.  This  noise  may  come  from  the  natural  environ¬ 
ment  of  the  system  or  can  be  added  to  the  system  artificially.  In  the  former 
case  (i.e.,  a  given  noise  structure),  determining  the  best  states  to  measure 
for  obtaining  the  best  precursor  is  an  important  problem  for  future  work. 
In  the  latter  setting,  we  can  consider  finding  the  optimal  location  for  in¬ 
jection  of  noise  to  result  in  the  clearest  possible  precursor. 

3.  In  some  systems,  the  naturally  occurring  perturbations  cannot  be  modeled 
as  white  noise,  but  rather  as  time  periodic  disturbances.  Precursors  and 
monitoring  systems  for  such  settings  deserve  investigation. 

4.  Experimental  testing  of  the  use  of  the  monitoring  systems  presented  here 
on  an  axial  flow  compression  system  rig  would  be  worthwhile. 
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Appendix  A 


Routh-Hurwitz  Calculation  for  Proposition  4.6 


Use  the  same  notation  as  in  Proposition  4.6  and  letting  7  =  a+ju ;  and  7  =  a—juj. 
Equation  (4.78)  becomes 

A6  —  2crA5  +  (2(c2  —  a)  +  a2  +  uj2)  A4  +  2a(2fl  —  e2)A3  + 

((c2  —  a)2  —  2a(a2  +  w2))A2  +  2(c2  —  a)aa  A  +  a2  (cr2  +  a;2)  =  0 


Applying  the  Routh-Hurwitz  criterion  to  the  equation  above,  we  obtain  the 
Routh  array 


1  2(c  —  a)  +  a 

—2a  2a  (2a  —  c ) 

c  +  a2  +  uj2  c2  —  ac  —  2  a(a2  +  uj2) 


(c  —  a)2  —  2a(a2  +  u>2)  a2(a2  +  uj2) 
2  aa(c-a)  0 


3  2  <Tc(a-(<r2+u2)) 

c+a'2+ci>2 

S2  A 

|  —2ac3a(c+a2+u2) 

S  a(a2+u2+c)-(a+c)2 

s°  a2(a2  +  uj2) 


2crac(c—  a+cr2+u2) 
c+a2+uj2 

a2(a2  +  uj2) 

0 

0 


a2(a2  +  uj2) 

0 

0 

0 

0 


0 

0 

0 

0 

0 
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where  A  =  ^2+^2)(^2+^(a+c)2(^2)  If  a  <  0  and  a  <  0.  then  all  entries  in 
the  first  column  of  this  array  are  positive.  Therefore,  under  this  condition  all 
solutions  of  (4.78)  have  negative  real  part. 
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Appendix  B 


Frequently  Used  Notation 


ck 


dx 
’  dt 


Rn 


C_ 


I 

eAt 

det(A) 

A 

Re(a ) 

7 

|ft| 

N{t),  n(t) 
<  X  > 


partial  derivatives  with  respect  to  x 

functions  k  times  differentiable 

time  derivative  of  x 

real  n  dimensional  space 

open  left  half  complex  plane 

parameters 

n  x  n  identity  matrix 

exponential  of  matrix  At 

determinant  of  n  x  n  matrix  A 

eigenvalue 

real  part  of  a 

complex  conjugate  of  7 

absolute  value  of  a 

white  Gaussian  noise 

expected  value  of  X 

time  average  over  all  time 
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■sgn(a)  sign  of  real  number  a 

A11  complex  conjugate  transpose  of  matrix  A 

g~1{x)  reciprocal  of  function  g(x) 
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